陸地棉自然群體纖維強(qiáng)度性狀的全基因組關(guān)聯(lián)分析

賈冰，馬建江，宋吉坤,，裴文鋒,,3，吳嫚，臧新山，張金發(fā)，陳全家*，于霽雯,,3*

（1.新疆農(nóng)業(yè)大學(xué)農(nóng)學(xué)院/ 棉花教育部工程研究中心，烏魯木齊830052；2.中國(guó)農(nóng)業(yè)科學(xué)院棉花研究所/ 棉花生物學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，河南 安陽455000；3.中國(guó)農(nóng)業(yè)科學(xué)院西部農(nóng)業(yè)研究中心，新疆 昌吉831100）

棉花作為世界上最重要的經(jīng)濟(jì)作物之一，可以提供大量的天然纖維[1]。 陸地棉因產(chǎn)量高、品質(zhì)較好在棉花生產(chǎn)中占主導(dǎo)地位[2]。 隨著人們生活水平的不斷提高， 對(duì)棉花纖維品質(zhì)的要求也越來越高。纖維強(qiáng)度是決定棉花纖維經(jīng)濟(jì)價(jià)值和紡織質(zhì)量的關(guān)鍵特性之一[3]。 纖維強(qiáng)度是數(shù)量性狀，由多基因共同調(diào)控，同時(shí)又受環(huán)境影響，很難通過傳統(tǒng)的育種技術(shù)加以改良[4]。 隨著分子生物學(xué)及生物信息技術(shù)的快速發(fā)展，分子標(biāo)記技術(shù)極大的推動(dòng)了數(shù)量性狀遺傳基礎(chǔ)的解析及改良。 因此，對(duì)棉花纖維強(qiáng)度與分子標(biāo)記進(jìn)行關(guān)聯(lián)分析以及鑒定纖維強(qiáng)度相關(guān)的優(yōu)異位點(diǎn)， 將進(jìn)一步加速棉纖維強(qiáng)度的改良進(jìn)程， 并為研究纖維強(qiáng)度形成的遺傳基礎(chǔ)提供理論依據(jù)。

棉花纖維的產(chǎn)量和質(zhì)量與棉花纖維發(fā)育進(jìn)程息息相關(guān)[5]。 棉纖維由胚珠表皮單個(gè)細(xì)胞分化而成，發(fā)育過程分為4 個(gè)時(shí)期：纖維起始分化、纖維伸長(zhǎng)期、次生壁加厚期以及成熟期[6]。纖維次生壁加厚期主要是纖維素的合成及沉淀的過程，纖維強(qiáng)度性狀主要受到該時(shí)期眾多基因的調(diào)控[7]。近年來，隨著分子生物技術(shù)的發(fā)展，特別是棉花全基因組測(cè)序的完成，關(guān)于棉花纖維強(qiáng)度關(guān)鍵位點(diǎn)、基因的挖掘及其功能的相關(guān)報(bào)道已逐漸增多。Deng 等[8]利用棉花多世代分離群體結(jié)合簡(jiǎn)單重復(fù)序列 （Simple sequence repeat，SSR）標(biāo)記，鑒定到17 個(gè)與纖維強(qiáng)度相關(guān)的數(shù)量性狀位點(diǎn) （Quantitative trait locus，QTL）。 Ma 等[9]利用陸地棉重組自交系群體及回交群體，鑒定50 個(gè)與纖維強(qiáng)度相關(guān)的QTL。Zhang 等[10]利用陸地棉重組自交系群體，構(gòu)建了全基因組高密度遺傳圖譜， 鑒定到33 個(gè)與纖維強(qiáng)度相關(guān)的穩(wěn)定的QTL。 Huang 等[11]利用CottonSNP 63K 芯片和503 份陸地棉種質(zhì)資源構(gòu)建高密度圖譜，鑒定到10個(gè)與纖維強(qiáng)度相關(guān)的QTL。 Ma 等[12]利用419 份陸地棉核心種質(zhì)，通過重測(cè)序技術(shù)進(jìn)行全基因組關(guān)聯(lián)分析， 檢測(cè)到533 個(gè)與纖維強(qiáng)度相關(guān)的顯著SNP位點(diǎn)，并鑒定了1 個(gè)與纖維強(qiáng)度相關(guān)的基因。 已有的研究獲得了大量的纖維強(qiáng)度相關(guān)的QTL 位點(diǎn)，但在多環(huán)境下穩(wěn)定的關(guān)鍵位點(diǎn)不多，仍需進(jìn)一步挖掘控制纖維強(qiáng)度性狀的關(guān)鍵位點(diǎn)和基因。

本研究以83 份纖維強(qiáng)度差異顯著的陸地棉為供試材料，利用63K 芯片對(duì)群體進(jìn)行基因型分析，得到15 369 個(gè)SNP（Single nuclotide polymorphism，單核苷酸多態(tài)性）位點(diǎn)，對(duì)5 個(gè)環(huán)境及最佳線性無偏估計(jì)值（Best linear unbiased prediction，BLUP）的纖維強(qiáng)度進(jìn)行全基因組關(guān)聯(lián)分析，解析棉花纖維強(qiáng)度性狀形成的遺傳基礎(chǔ)。

1 材料與方法

1.1 試驗(yàn)材料及種植方法

2014―2016 年分別在河南省安陽市安陽縣（14_Ay、15_Ay、16_Ay）、2016 年 新 疆 維 吾 爾 自 治區(qū)阿拉爾市十團(tuán) （16_Ale）、2016 年海南省三亞市吉陽區(qū)（16_Sy）等5 個(gè)環(huán)境種植83 份陸地棉材料，材料信息及種植模式參考文獻(xiàn)13、14。

1.2 表型測(cè)定與分析

在棉花成熟時(shí)， 選擇棉株中部相同部位的20個(gè)棉鈴，統(tǒng)一用皮輥軋花機(jī)軋花（SY-50）。纖維斷裂比強(qiáng)度由農(nóng)業(yè)部棉花品質(zhì)檢驗(yàn)監(jiān)督檢測(cè)中心的HFT9000 儀器檢測(cè)完成。

利用SPSS 25.0 對(duì)83 份陸地棉材料在不同環(huán)境的纖維強(qiáng)度進(jìn)行描述性統(tǒng)計(jì)分析和相關(guān)性分析，使用QTL IciMapping 4.2.0 進(jìn)行廣義遺傳力分析，使用R 語言計(jì)算BLUP 值。

1.3 纖維強(qiáng)度性狀QTL 定位

利用課題之前使用CottonSNP 63K 芯片對(duì)83份陸地棉材料進(jìn)行基因型分析得到的15 369 個(gè)SNP 位點(diǎn)[13]，使用TASSEL 5.0 中的廣義線性模型（General linear model,GLM） 對(duì)纖維強(qiáng)度性狀進(jìn)行全基因組關(guān)聯(lián)分析。 GWAS 顯著位點(diǎn)的閾值（Bonferroni 法）P為6.51E-05[15-17]，即將－lg（P）≥4.19 的位點(diǎn)作為與纖維強(qiáng)度顯著關(guān)聯(lián)的SNP 位點(diǎn)，使用R 語言繪制曼哈頓圖。 本試驗(yàn)認(rèn)為：在3個(gè)及以上不同環(huán)境檢測(cè)到的同一位置的SNP 是穩(wěn)定的SNP 位點(diǎn)。 根據(jù)每條染色體特定的LD（Linkage disequilibrium，連鎖不平衡）的衰減距離[13]，將穩(wěn)定的SNP 位點(diǎn)物理位置±LD 值作為QTL 區(qū)間，重疊的區(qū)間合并為1 個(gè)QTL 區(qū)間，QTL 命名方式為：q-性狀-染色體-QTL 編號(hào)[18]。

1.4 候選基因分析

將QTL 區(qū)間內(nèi)的基因作為候選基因庫， 根據(jù)海島棉海7124 和陸地棉TM-1 的轉(zhuǎn)錄組數(shù)據(jù)庫[19]（https://www.ncbi.nlm.nih.gov/bioproject/PRJNA490626/）以及擬南芥同源基因功能注釋，篩選與纖維強(qiáng)度相關(guān)基因作為候選基因。

2 結(jié)果與分析

2.1 陸地棉纖維強(qiáng)度表型分析

通過對(duì)5 個(gè)環(huán)境83 份陸地棉材料的纖維斷裂比強(qiáng)度進(jìn)行統(tǒng)計(jì)分析，發(fā)現(xiàn)纖維斷裂比強(qiáng)度性狀具有廣泛的表型變異。 該群體纖維斷裂比強(qiáng)度為22.50～38.50 cN·tex－1，在14_Ay、15_Ay、16_Ay 三個(gè)環(huán)境中纖維斷裂比強(qiáng)度的平均值分別為28.94 cN·tex－1、28.70 cN·tex－1和29.58 cN·tex－1， 而在16_Sy 和16_Ale 兩個(gè)環(huán)境中纖維斷裂比強(qiáng)度的平均值分別為26.34 cN·tex－1和25.47 cN·tex－1（表1和圖1）。 統(tǒng)計(jì)結(jié)果表明，在安陽種植的棉花材料的纖維斷裂比強(qiáng)度3 年的分布趨于一致，與16_Sy 和16_Ale 呈現(xiàn)一定的差異。 在5 個(gè)環(huán)境中纖維斷裂比強(qiáng)度的變異系數(shù)為5.55%～8.44%， 廣義遺傳力高達(dá)88.67%， 說明該棉花群體纖維斷裂比強(qiáng)度變異較豐富，遺傳相對(duì)穩(wěn)定。 相關(guān)性分析（表2）表明，纖維斷裂比強(qiáng)度在5 個(gè)環(huán)境中相關(guān)系數(shù)為0.6～0.76，均呈顯著正相關(guān)。 對(duì)5 個(gè)環(huán)境的纖維斷裂比強(qiáng)度表型數(shù)據(jù)進(jìn)行正態(tài)分布檢驗(yàn)， 發(fā)現(xiàn)除16_Ale外，其他環(huán)境的纖維斷裂比強(qiáng)度性狀基本符合正態(tài)分布，表明該群體適合進(jìn)行棉花纖維斷裂比強(qiáng)度的全基因組關(guān)聯(lián)分析。 同時(shí)，83 份材料在安陽種植纖維斷裂比強(qiáng)度較大，而在新疆和海南種植則較小，說明纖維斷裂比強(qiáng)度受到環(huán)境的影響， 需要對(duì)單個(gè)環(huán)境進(jìn)行全基因組關(guān)聯(lián)分析。

2.2 全基因組關(guān)聯(lián)分析及纖維斷裂比強(qiáng)度QTL

利用TASSEL 軟件的GLM（Q）模型對(duì)14_Ay、15_Ay、16_Ay、16_Sy、16_Ale 5 個(gè)環(huán)境的纖維斷裂比強(qiáng)度及BLUP 進(jìn)行全基因組關(guān)聯(lián)分析（圖2），分別得到39、134、14、7、16 和55 個(gè)顯著相關(guān)的SNP位點(diǎn)，共計(jì)265 個(gè)SNP 位點(diǎn)。 去冗余后得到190 個(gè)SNP 位點(diǎn)，分布于除A13 外的其他所有染色體上，

對(duì)纖維斷裂比強(qiáng)度的表型解釋為18%～40%，其中121 個(gè)顯著SNP 位點(diǎn)位于Dt 亞基因組，69 個(gè)顯著SNP 位點(diǎn)位于At 亞基因組。 在3 個(gè)及以上環(huán)境中出現(xiàn)的穩(wěn)定的SNP 位點(diǎn)有19 個(gè)，結(jié)合其染色體的LD 的衰減距離[13]，得到9 個(gè)QTL 區(qū)間，分別位于A01、A06、D05、D08、D10、D11 和D13 等7 條染色體上 （表3）， 其中4 個(gè)QTL 區(qū)間與前人定位的QTL 位置重疊，其它5 個(gè)QTL 區(qū)間未見報(bào)道。

表1 陸地棉群體纖維斷裂比強(qiáng)度表型統(tǒng)計(jì)

表2 五個(gè)環(huán)境群體纖維斷裂比強(qiáng)度的相關(guān)性分析

圖1 五個(gè)環(huán)境下群體纖維斷裂比強(qiáng)度的頻率分布直方圖

2.3 纖維斷裂比強(qiáng)度候選基因挖掘

纖維斷裂比強(qiáng)度主要在纖維次生壁加厚期受纖維素合成的影響， 且海島棉海7124 和陸地棉TM-1 的纖維斷裂比強(qiáng)度存在顯著差異， 所以篩選在20～25 DPA纖維中優(yōu)勢(shì)表達(dá)， 且在海7124 和TM-1 間差異表達(dá)，基因功能注釋與纖維斷裂比強(qiáng)度相關(guān)的基因作為候選基因。 9 個(gè)穩(wěn)定QTL 區(qū)間包含1 177 個(gè)基因，結(jié)合海7124 和TM-1 轉(zhuǎn)錄組數(shù)據(jù)，鑒定到4 個(gè)與細(xì)胞壁形成或纖維素合成相關(guān)的基因（表4），推測(cè)這些基因可能影響纖維斷裂比強(qiáng)度。Gh_D13G2032編碼產(chǎn)物為果膠甲酯酶， 在海7124 開花后25 d 的纖維中表達(dá)量最高， 顯著高于其在TM-1 同時(shí)期纖維中的表達(dá)量。Gh_A01G1734在海7124 開花后20 d 的纖維中表達(dá)量最高，顯著高于其在TM-1 中的表達(dá)量。Gh_D10G0333編碼纖維素合酶，Gh_D11G1626屬于COBL基因家族，這2 個(gè)基因都是在20 DPA 的纖維中表達(dá)量最高，且在海7124 的表達(dá)量明顯高于TM-1（圖3）。

圖2 陸地棉纖維斷裂比強(qiáng)度全基因組關(guān)聯(lián)分析

表3 纖維斷裂比強(qiáng)度相關(guān)QTL 信息

表3 （續(xù)）

表4 纖維斷裂比強(qiáng)度相關(guān)候選基因信息

圖3 海7124 和TM-1 纖維發(fā)育不同時(shí)期4 個(gè)候選基因的表達(dá)模式

3 討論

棉花纖維品質(zhì)主要看“長(zhǎng)”、“強(qiáng)”、“細(xì)” ，其中“強(qiáng)”就是纖維斷裂比強(qiáng)度，對(duì)棉纖維的加工和產(chǎn)品的質(zhì)量有直接的影響。 因此，棉花纖維斷裂比強(qiáng)度改良一直是棉花纖維品質(zhì)遺傳改良的熱點(diǎn)。

本研究中用于全基因組關(guān)聯(lián)分析的83 份陸地棉材料，具有豐富的表型變異。 通過對(duì)5 個(gè)環(huán)境的纖維斷裂比強(qiáng)度表型數(shù)據(jù)及BLUP 值進(jìn)行全基因組關(guān)聯(lián)分析，共得到了265 個(gè)顯著SNP 位點(diǎn)，非冗余的顯著位點(diǎn)190 個(gè)， 其中位于Dt 亞基因組上的SNP 位點(diǎn)是At 亞基因組的顯著SNP 位點(diǎn)的1.7倍。 前人研究表明Dt 亞基因組對(duì)纖維品質(zhì)性狀的遺傳貢獻(xiàn)比At 亞基因組更大， 本研究結(jié)果與前人的研究結(jié)果一致[19]。 另外，19 個(gè)SNP 位點(diǎn)出現(xiàn)在3個(gè)及以上環(huán)境中， 合并為9 個(gè)QTL 區(qū)間， 其中qFS-A01-1包含Sun 等[21]檢測(cè)到的顯著SNP 位點(diǎn)，qFS-D05-1與qFS-D13-2均與Wang 等[22]檢測(cè)到的QTL 區(qū)間重疊，qFS-D10-1與Gu 等[23]檢測(cè)到的區(qū)間重疊。 此外，本研究還發(fā)現(xiàn)了5 個(gè)新的穩(wěn)定QTL區(qū)間，為棉花纖維斷裂比強(qiáng)度及分子標(biāo)記輔助選擇育種提供參考。

在9 個(gè)QTL 區(qū)間內(nèi)篩選了22 個(gè)在纖維發(fā)育次生壁加厚期優(yōu)勢(shì)表達(dá)、 且在海島棉與陸地棉之間差異表達(dá)的基因， 并鑒定到4 個(gè)可能影響纖維斷裂比強(qiáng)度的候選基因。 其中Gh_D13G2032的編碼產(chǎn)物為果膠甲酯酶， 催化細(xì)胞壁中果膠主鏈發(fā)生脫甲酯化反應(yīng)， 對(duì)于調(diào)節(jié)植物中的細(xì)胞膨脹非常重要。 在前人研究中，海島棉纖維比陸地棉更細(xì)更長(zhǎng)可能是受果膠甲酯酶含量差異的影響[24]。Gh_D11G1626作為COBL基因家族成員，編碼植物特異性糖基磷脂酰肌醇錨定蛋白， 該蛋白通過影響微纖絲在細(xì)胞壁的正確排列而參與植物纖維素的合成，在纖維發(fā)育過程中影響纖維素的沉淀，進(jìn)而影響棉花纖維斷裂比強(qiáng)度[25]。Gh_A01G1734的同源基因AT5G06270顯著影響擬南芥的根毛數(shù)量， 且與擬南芥的表皮毛發(fā)育密切相關(guān)[26]。Gh_D10G0333的功能注釋為纖維素合酶GhCESA8， 在棉花纖維發(fā)育過程中纖維素的快速大量積累起到了增強(qiáng)作用[27]。

4 結(jié)論

本研究對(duì)83 份棉花材料在5 個(gè)環(huán)境的纖維斷裂比強(qiáng)度及其BLUP 值進(jìn)行全基因組關(guān)聯(lián)分析，共檢測(cè)到19 個(gè)與纖維斷裂比強(qiáng)度有關(guān)的穩(wěn)定的SNP位點(diǎn)，得到9 個(gè)QTL 區(qū)間，其中5 個(gè)是本研究新發(fā)現(xiàn)的。并鑒定到4 個(gè)與纖維斷裂比強(qiáng)度有關(guān)的候選基因，為進(jìn)一步挖掘與棉花纖維斷裂比強(qiáng)度性狀相關(guān)的基因提供了理論依據(jù)。