整合GWAS和WGCNA篩選鑒定甘藍(lán)型油菜生物產(chǎn)量候選基因

王艷花 劉景森 李加納

整合GWAS和WGCNA篩選鑒定甘藍(lán)型油菜生物產(chǎn)量候選基因

王艷花1,2,**劉景森1,2,**李加納1,2,*

1西南大學(xué)農(nóng)學(xué)與生物科技學(xué)院/ 油菜工程研究中心, 重慶 400715;2西南大學(xué)現(xiàn)代農(nóng)業(yè)科學(xué)研究院, 重慶 400715

生物產(chǎn)量是作物獲得高產(chǎn)的重要基礎(chǔ), 對(duì)于甘藍(lán)型油菜(L.)尤其重要。本研究利用588份甘藍(lán)型油菜材料構(gòu)成的自然群體2年生物產(chǎn)量表型數(shù)據(jù)的全基因組關(guān)聯(lián)分析, 再結(jié)合高生物產(chǎn)量材料‘CQ45’和低生物產(chǎn)量材料‘CQ46’的轉(zhuǎn)錄組測序(RNA-seq)結(jié)果, 整合了6個(gè)甘藍(lán)型油菜材料6個(gè)部位(莖稈、葉片、花后30 d主軸與側(cè)枝種子、花后30 d主軸與側(cè)枝角果皮)的轉(zhuǎn)錄組數(shù)據(jù)構(gòu)建的加權(quán)共表達(dá)網(wǎng)絡(luò)分析(WGCNA), 篩選出與生物產(chǎn)量相關(guān)的候選基因。通過相關(guān)分析發(fā)現(xiàn), 2年間甘藍(lán)型油菜自然群體中生物產(chǎn)量對(duì)大多數(shù)產(chǎn)量相關(guān)性狀都具有正向效應(yīng); 自然群體2年生物產(chǎn)量分析的最佳模型均為K+PCA模型, 共檢測到9個(gè)顯著位點(diǎn)(< 1/385691或< 0.05/385691); 根據(jù)CQ45和CQ46共36組轉(zhuǎn)錄組數(shù)據(jù), 選擇MAD值為前5%的基因共計(jì)5052個(gè)用于構(gòu)建WGCNA, 通過篩選合并共得到了15個(gè)模塊, 其中5個(gè)基因共表達(dá)模塊分別與葉片、莖稈和花后30 d種子顯著性相關(guān); 整合了WGCNA中關(guān)鍵模塊的hub gene、GWAS分析得到的顯著SNP位點(diǎn)和極端表型差異基因確定候選基因, 它們的擬南芥同源基因?yàn)?、、? 這些基因在光合作用的卡爾文循環(huán)、碳同化、物質(zhì)積累等方面發(fā)揮重要作用。

GWAS; WGCNA; RNA Seq; 甘藍(lán)型油菜

油菜作為四大油料作物(大豆、油菜、向日葵、花生)之一, 由于其較強(qiáng)的適應(yīng)性, 在世界各地廣泛種植[1-3]。隨著人類對(duì)油菜的需求不斷增加, 油菜總產(chǎn)量增速不斷加快[4], 世界油料作物中油菜的種植面積及總產(chǎn)量已位居第二, 僅次于大豆[5-6]。同時(shí)油菜也是我國主要農(nóng)作物, 在作物種植面積中排名第五, 僅次于水稻、玉米、小麥、大豆[7-12]。近年來, 油菜多功能開發(fā)利用的持續(xù)發(fā)展, 形成了油菜的油用、菜用、花用、蜜用、飼用等“一菜多用”的模式, 大幅度的提高了油菜的種植效益[13]。在實(shí)際生產(chǎn)中, 我國油菜產(chǎn)量仍較低, 生產(chǎn)成本卻很高, 國際競爭力嚴(yán)重不足[14-15]。目前, 我國油菜超過60%依賴進(jìn)口, 并且種植面積和單位面積產(chǎn)量都停滯不前, 這對(duì)我國糧食安全戰(zhàn)略有著嚴(yán)重威脅[16]。因此, 提高油菜產(chǎn)量是我國油菜育種不斷追求的目標(biāo)[17-24]。

生物產(chǎn)量指單位面積土地上獲得的不包括根系的作物干物質(zhì)的總量, 它是描述莖、葉制造并積累有機(jī)物能力強(qiáng)弱的量[25]。生物產(chǎn)量體現(xiàn)作物總的生產(chǎn)力, 或者說作物光合作用碳同化的能力。多項(xiàng)報(bào)道已經(jīng)表明, 生物產(chǎn)量與經(jīng)濟(jì)產(chǎn)量密切相關(guān), 作物提高產(chǎn)量的方法通常有2種: 即在生物產(chǎn)量不變的情況下, 提高收獲指數(shù); 或收獲指數(shù)不變的情況下, 提高生物產(chǎn)量[26-30]。近年來, 高通量測序技術(shù)的不斷突破、成本不斷降低, 多樣本的轉(zhuǎn)錄組測序逐漸被應(yīng)用于系統(tǒng)研究生命科學(xué)問題[31], 傳統(tǒng)的少量樣本比較分析已經(jīng)無法有效處理海量的生物信息數(shù)據(jù)。因此, 生物信息分析工具和方法應(yīng)運(yùn)而生[32]。其中, 加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(weighted gene co-expression network analysis, WGCNA)能夠特異篩選出與目標(biāo)性狀高度關(guān)聯(lián)的基因, 并進(jìn)行模塊(module)化分類, 得到具有高度生物學(xué)意義的共表達(dá)模塊(co-expression module), 能夠有效篩選到核心基因(hub genes)[33]。WGCNA算法作為一種精準(zhǔn)、高效的生物信息學(xué)及生物數(shù)據(jù)挖掘方法, 被廣泛應(yīng)用于生物學(xué)各個(gè)領(lǐng)域, 利用WGCNA對(duì)馬鈴薯C16和C119品種構(gòu)建抗逆生理性狀高度關(guān)聯(lián)的權(quán)重網(wǎng)絡(luò), 為深入探究馬鈴薯抗逆分子機(jī)制提供新的數(shù)據(jù)資源[34-35]。本研究以生物產(chǎn)量性狀為主要研究對(duì)象, 通過GWAS分析, 鑒定出與生物產(chǎn)量性狀緊密相關(guān)的SNP位點(diǎn); 然后提取SNP位點(diǎn)前后共500 kb的序列片段作為候選區(qū)段, 提取候選區(qū)段內(nèi)所有基因并結(jié)合轉(zhuǎn)錄組分析結(jié)果篩選獲得差異候選基因; 同時(shí), 結(jié)合轉(zhuǎn)錄組數(shù)據(jù)構(gòu)建WGCNA共表達(dá)網(wǎng)絡(luò), 篩選鑒定出甘藍(lán)型油菜生物產(chǎn)量關(guān)鍵候選基因。

1 材料與方法

1.1 供試材料

以來自世界各地的588份油菜自交系作為材料, 其中455份來自中國(大部分來自重慶、湖南、湖北等地), 11份來自亞洲其他地區(qū), 102份來自歐洲, 13份來自北美洲, 7份來自澳大利亞。所有材料均由重慶市油菜工程技術(shù)研究中心提供。

1.2 材料種植及農(nóng)藝性狀考察

588份材料于2016、2018年種植于重慶北碚歇馬油菜基地(29o45′39.99″, 106o22′38.47″, 海拔238.57 m)。每年播種時(shí)間為9月20日左右, 按完全隨機(jī)區(qū)組設(shè)計(jì), 2個(gè)重復(fù)同時(shí)種植, 每個(gè)材料種植2行, 每行15株, 行距40 cm, 株距20 cm。田間管理選擇常規(guī)生產(chǎn)方式, 2個(gè)重復(fù)按照相同的管理措施和施肥水平進(jìn)行, 收獲時(shí)間為次年5月5日左右。于成熟期, 在每個(gè)小區(qū)中選擇5株長勢(shì)一致的植株。將植株分成上部分枝和下部莖稈兩部分, 并將上部分枝裝進(jìn)網(wǎng)袋中。上部分枝和下部莖稈分別徹底曬干后, 使用天平秤進(jìn)行稱量, 將兩部分重量相加得到植株的總干重。每份材料的5個(gè)樣本干重的平均值為該材料的生物產(chǎn)量(biomass yield, BY)。

1.3 表型數(shù)據(jù)分析

利用Microsoft Excel軟件初步整理588份自然群體的表型數(shù)據(jù), 并計(jì)算單個(gè)環(huán)境及差值條下生物產(chǎn)量的平均值、標(biāo)準(zhǔn)差、變異系數(shù)等, 并用SPSS繪制生物產(chǎn)量的正態(tài)分布圖; 利用SPSS統(tǒng)計(jì)分析軟件對(duì)自然群體進(jìn)行各性狀進(jìn)行相關(guān)分析; 利用R軟件對(duì)多年環(huán)境下自然群體的生物產(chǎn)量表型數(shù)據(jù)進(jìn)行最佳線性無偏預(yù)測(best linear unbiased prediction, BLUP)處理。

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

利用本實(shí)驗(yàn)室已有588份甘藍(lán)型油菜重測序的數(shù)據(jù), 通過基因型分析, 最終獲得385,691個(gè)可利用的SNP標(biāo)記數(shù)據(jù)[36]。本研究利用這些標(biāo)記結(jié)合重慶2017、2019年和BLUP環(huán)境的生物產(chǎn)量(BY)表型數(shù)據(jù)行全基因組關(guān)聯(lián)分析。本研究使用6種統(tǒng)計(jì)模型在Tassel 5.2.1[37]軟件中執(zhí)行對(duì)甘藍(lán)型油菜生物產(chǎn)量的關(guān)聯(lián)分析, 即: Q模型、naive模型、K模型、PCA模型、Q+K模型和PCA+K模型。Q代表群體結(jié)構(gòu), PCA代表主成分, K代表親緣關(guān)系; naive、Q和PCA模型在一般線性模型(general linear model, GLM)下運(yùn)行, K、Q+K和PCA+K模型在混合線性模型(mixed linear model, MLM)下運(yùn)行。利用SAS軟件對(duì)這6種模型的運(yùn)算結(jié)果中的-lg()的觀測值和-lg()期望值繪制QQ (Quantile-quantile)圖。通過比較QQ圖確定最佳模型, 并選擇最佳模型下生物產(chǎn)量性狀的全基因組關(guān)聯(lián)分析結(jié)果作為進(jìn)一步分析對(duì)象。本試驗(yàn)使用的SNP數(shù)據(jù)為385,691個(gè), 因此規(guī)定值小于閾值(1/385,691=2.593E-06和0.05/385,691=1.296E-07)的位點(diǎn)為顯著關(guān)聯(lián)位點(diǎn)。利用Haploview[38]軟件制作Manhattan圖, 將通過關(guān)聯(lián)分析檢測到的染色體組上與各性狀顯著相關(guān)的標(biāo)記位點(diǎn)可視化。

1.5 轉(zhuǎn)錄組測序分析

從供試群體中選擇1個(gè)高生物產(chǎn)量(CQ45)和1個(gè)低生物產(chǎn)量(CQ46)材料, 分別選取蕾苔期莖稈(J)和葉片(Le)、采用掛花標(biāo)記的方法選取開花后30 d主軸種子(30S)、主軸角果皮(30P)、側(cè)枝種子(30CS)和側(cè)枝角果皮(30CP) 6個(gè)時(shí)期的樣品, 每個(gè)樣品取2個(gè)生物重復(fù), 送北京諾禾致源生物信息科技有限公司進(jìn)行轉(zhuǎn)錄組測序, 在OmicShare (基迪奧)云平臺(tái)完成(http://www.omicshare.com) GO和KEGG富集分析, 篩選差異表達(dá)基因。

1.6 加權(quán)共表達(dá)網(wǎng)絡(luò)分析

使用R軟件包preprocessCore中的normalize. quantiles函數(shù)對(duì)轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化。然后計(jì)算每個(gè)基因的MAD (median absolute deviation), 選擇的MAD值為前5%的基因作為樣本間變異大的基因群用于WGCNA構(gòu)建。使用average-linkage層次聚類法對(duì)基因進(jìn)行聚類, 按照混合動(dòng)態(tài)剪切樹的標(biāo)準(zhǔn), 并設(shè)置每個(gè)基因網(wǎng)絡(luò)模塊最少的基因數(shù)為30。在使用動(dòng)態(tài)剪切法在確定基因模塊后, 我們依次計(jì)算每個(gè)模塊的特征向量值(eigengenes), 然后對(duì)模塊進(jìn)行聚類分析, 將距離較近的模塊合并成新的模塊, 設(shè)置height = 0.25、deepSplit = 3、minModuleSize = 30。計(jì)算所得到模塊的特征向量與樣本來源組織的相關(guān)性, 選擇相關(guān)系數(shù)>0.6的模塊作為與組織部位顯著相關(guān)模塊。通過特征向量基因分析確定目標(biāo)基因模塊, 對(duì)所有的模塊中具有代表性的基因-特征向量基因(module eigengene, ME)進(jìn)行聚類分析, 進(jìn)一步進(jìn)行兩兩模塊之間的ME相關(guān)性分析, ME之間的相關(guān)性越高, 其所在的模塊的相關(guān)性也就越高; 將模塊中的所有基因和ME基因在所有樣本中分別進(jìn)行表達(dá)水平分析, GO和KEGG分析。最后使用Cytoscape[39]將所關(guān)注模塊的ME基因的共表達(dá)網(wǎng)絡(luò)進(jìn)行可視化。

2 結(jié)果與分析

2.1 表型數(shù)據(jù)統(tǒng)計(jì)及各性狀的相關(guān)性分析

對(duì)2017年和2019年588份自然群體的生物產(chǎn)量的表型變異情況進(jìn)行統(tǒng)計(jì)分析發(fā)現(xiàn), 在不同環(huán)境下自群體的生物產(chǎn)量呈連續(xù)變異, 近似正態(tài)分布, 說明生物產(chǎn)量是受多基因控制的數(shù)量性狀, 符合GWAS分析的要求。各產(chǎn)量相關(guān)性狀的相關(guān)分析結(jié)果表明, 每角果粒數(shù)在2017年和2019年與莖稈干重、經(jīng)濟(jì)產(chǎn)量、收獲指數(shù)和生物產(chǎn)量均呈正相關(guān), 而莖稈干重僅在2017年呈現(xiàn)出顯著正相關(guān); 莖稈干重與上部生物產(chǎn)量、經(jīng)濟(jì)產(chǎn)量、生物產(chǎn)量在2年間均呈現(xiàn)顯著正相關(guān); 上部生物產(chǎn)量與經(jīng)濟(jì)產(chǎn)量、收獲指數(shù)、生物產(chǎn)量在2年間均呈現(xiàn)顯著正相關(guān); 而經(jīng)濟(jì)產(chǎn)量和收獲指數(shù)與生物產(chǎn)量在2年間同樣呈現(xiàn)正相關(guān)。多數(shù)性狀間的相關(guān)性均達(dá)到顯著或極顯著水平, 表明產(chǎn)量相關(guān)的各性狀間存在顯著相關(guān)作用(表1)。

2.2 甘藍(lán)型油菜生物產(chǎn)量性狀的GWAS分析

各生物產(chǎn)量性狀在6種模型下的關(guān)聯(lián)分析后得到的結(jié)果進(jìn)行Quantile-Quantile散點(diǎn)圖(QQ Plot)的繪制, 選擇與預(yù)測值最接近的一種模型作為最佳模型, 2017年、2019年和BLUP環(huán)境下最佳模型為K+PCA模型(圖1)。在最佳模型下, 利用本實(shí)驗(yàn)室已有的385,691個(gè)SNP, 以值小于閾值(1/385,692 = 2.593E-06)確定顯著關(guān)聯(lián)SNP位點(diǎn), 并繪制曼哈頓圖(圖2)。生物產(chǎn)量2017、2019年和BLUP環(huán)境下共檢測到9個(gè)顯著關(guān)聯(lián)標(biāo)記位點(diǎn)(表2)。其中, 在2017年環(huán)境中檢測到7個(gè)顯著關(guān)聯(lián)標(biāo)記位點(diǎn), 分別位于A03、A07、C01 (2個(gè))、C04 (2個(gè))、C06染色體, 并且位于同一條染色體的顯著關(guān)聯(lián)位點(diǎn)其置信區(qū)間高度重合; 2019年檢測到2個(gè)顯著關(guān)聯(lián)位點(diǎn), 分別位于A07和C09染色體; 而BLUP環(huán)境下共檢測到4個(gè)顯著關(guān)聯(lián)位點(diǎn), 分別位于A07、C03、C06 (2個(gè)), 位于C06染色體的2個(gè)顯著關(guān)聯(lián)位點(diǎn)其置信區(qū)間高度重合。9個(gè)SNP的貢獻(xiàn)率在5.64%~7.98%, 其中2017年檢測到的位于C06染色體的SNP (S16_ 26169577)貢獻(xiàn)率最高。

表1 自然群體2年各產(chǎn)量相關(guān)性狀間的相關(guān)分析

*、**分別表示在0.05和0.01水平相關(guān)性顯著(雙尾)。.

*,**: significant correlation at the 0.05 and 0.01 probability levels, respectively. SWSI: seeds weight per silique index; TSW: 1000-grain weight; SNPS: seeds number per silique; ST: stem dry weight; CBY: canopy biomass yield; SY: seed yield; HI: harvest index; BY: biomass yield.

2.3 極端表型材料RNA-Seq分析

在葉片中檢測到極端材料的差異基因共6820個(gè), CQ45對(duì)CQ46顯著上調(diào)表達(dá)的基因有3516個(gè), 顯著下調(diào)表達(dá)的基因有3304個(gè)(圖3-A)。生物分子功能類中相關(guān)性最高的是氧化還原酶活性(GO:0016903)、催化活性(GO:0003824)等, 而參與基因數(shù)量較多的生物過程是單體代謝過程(GO:0044710)、光合作用(GO:0015979)等(圖3-B); 葉片中差異基因主要集中在碳水化合物代謝途徑(圖3-C)。

花后30 d主軸角果皮中檢測到17,309個(gè)顯著差異表達(dá)基因, CQ45對(duì)于CQ46有2682個(gè)顯著上調(diào)表達(dá)基因, 3570個(gè)顯著下調(diào)基因(圖3-A)。差異基因中分子功能相關(guān)性較高的是葡萄酶轉(zhuǎn)移(GO:0046527)等, 生物過程相關(guān)性最高的有應(yīng)激反應(yīng)(GO:0006950)、單體生物合成過程(GO:0044711)等(圖3-B); 相關(guān)性較高的途徑是氨基酸的生物合成以及脂質(zhì)代謝途徑(圖3-C)。

花后30 d側(cè)枝角果皮中檢測到5431個(gè)顯著差異表達(dá)基因, 其中CQ45對(duì)于CQ46有2173個(gè)顯著上調(diào)表達(dá)基因, 3258個(gè)顯著下調(diào)基因(圖3-A)。側(cè)枝角果皮中分子功能相關(guān)性最高的是營養(yǎng)庫活性(GO:0045735)等, 參與的生物過程中相關(guān)性較高的是對(duì)含氧化合物的反應(yīng)(GO:1901700)等(圖3-B); 這些差異基因相關(guān)性最高的是能量代謝途徑(圖3-C)。

花后30 d主軸種子中檢測到6948個(gè)顯著差異表達(dá)基因, CQ45對(duì)于CQ46有4080個(gè)顯著上調(diào)表達(dá)基因, 2868個(gè)顯著下調(diào)基因(圖3-A)。分子功能相關(guān)性最高的是微管運(yùn)動(dòng)活性(GO:0003777)、硫葡萄糖苷酶活性(GO:0019137)等, 生物過程相關(guān)性較高的是細(xì)胞循環(huán)(GO:0007049)、細(xì)胞分裂(GO:0051301) (圖3-B); 相關(guān)性較高的代謝途徑是碳水化合物代謝中的淀粉和蔗糖代謝(圖3-C)。

花后30 d側(cè)枝種子中檢測到17,309個(gè)顯著差異表達(dá)基因, CQ45對(duì)于CQ46有11489個(gè)顯著上調(diào)表達(dá)基因, 5820個(gè)顯著下調(diào)基因(圖3-A)。分子功能相關(guān)性最高的是微管運(yùn)動(dòng)活性(GO:0003777)、營養(yǎng)庫活性(GO:0045735)硫葡萄糖苷酶活性(GO:0019137)等, 參與的生物過程相關(guān)性較高的是細(xì)胞循環(huán)(GO:0007049)等, 其結(jié)果與30 d主軸角果種子比較相似(圖3-B); 相關(guān)性較高的是氨基酸代謝中的丙氨酸、天冬氨酸和谷氨酸代謝(圖3-C)。

表2 最佳模型下生物產(chǎn)量顯著位點(diǎn)表

(圖3)

(圖3)

(圖3)

A: 極端表型材料各組織差異基因數(shù)量; B: 極端表型材料各組織差異基因的GO富集分析; C: 極端表型材料各組織差異基因的KEGG分析。

A: the number of differentially expressed genes in tissues of extreme phenotypic materials; B: GO enrichment analysis of differentially expressed genes in tissues of extreme phenotypic materials; C: KEGG analysis of differentially expressed genes in tissues of extreme phenotype materials.

莖稈中檢測到極端表型材料的差異基因12,867個(gè), CQ45對(duì)于CQ46有6452個(gè)顯著上調(diào)表達(dá)基因, 6415個(gè)顯著下調(diào)基因(圖3-A)。分子功能相關(guān)性最高的是轉(zhuǎn)移酶活性(GO:0016758)、UDP葡萄糖基轉(zhuǎn)移酶活性(GO:0035251)等, 參與較多的生物過程是植物型次生細(xì)胞壁生物發(fā)生(GO:0009834)、細(xì)胞壁合成(GO:0042546)等(圖3-B); 參與相關(guān)性較高的代謝途徑是到次生代謝產(chǎn)物的生物合成等途徑(圖3-C)。

2.4 加權(quán)共表達(dá)網(wǎng)絡(luò)(WGCNA)分析

2.4.1 基因共表達(dá)網(wǎng)絡(luò)的模塊生成 根據(jù)36個(gè)樣本得到的RNA-Seq的結(jié)果, 計(jì)算各基因的MAD值, 取MAD為前5%共計(jì)5052個(gè)基因, 進(jìn)行基因表達(dá)水平的層次聚類(圖4-A); 使用R軟件包WGCNA構(gòu)建權(quán)重共表達(dá)網(wǎng)絡(luò), 選取擬合曲線第1次接近0.9時(shí)的閾值參數(shù)β (β=10)篩選共表達(dá)模塊(圖4-B, C),共表達(dá)網(wǎng)絡(luò)符合無尺度網(wǎng)絡(luò)。使用動(dòng)態(tài)剪切法在確定基因模塊后, 依次計(jì)算每個(gè)模塊的特征向量值, 然后對(duì)模塊進(jìn)行聚類分析, 將距離較近的模塊合并成新的模塊后共得到了15個(gè)模塊(圖4-D), 聚類樹中每個(gè)枝代表1個(gè)模塊, 分別是Black、Cyan、Green、Magenta、Pink、Purple、Red、Salmon、Tan、Turquoise、Yellow、Greenyellow、Brown、Blue模塊, Grey模塊是無法聚集到其他模塊的基因集合(包含11個(gè)基因)。在15個(gè)模塊中基因數(shù)量33~1843不等, 其中包含基因數(shù)量最多的是Turquoise模塊(1843個(gè)), 最少的是Cyan模塊(33個(gè)), 各個(gè)基因之間的表達(dá)關(guān)系如圖4-E。

2.4.2 模塊與組織之間相關(guān)性 本研究分別計(jì)算這15個(gè)模塊的特征向量與樣本來源組織的相關(guān)性, 得到5個(gè)顯著相關(guān)的模塊, 其中與葉片(L)極顯著相關(guān)的模塊是Turquoise模塊, 相關(guān)系數(shù)達(dá)到0.77; 與莖稈顯著相關(guān)的模塊為Green、Magenta、Pink模塊, 其相關(guān)系數(shù)分別為0.9、0.62、0.75; 與30 d主軸種子顯著相關(guān)的模塊是Yellow模塊, 其相關(guān)系數(shù)為0.6 (圖5-A)。

通過特征向量分析確定目標(biāo)基因模塊, 對(duì)所有的模塊中具有代表性的基因-特征向量基因(module eigengene, ME)進(jìn)行聚類分析(圖5-B), ME之間的相關(guān)性越高, 其所在的模塊的相關(guān)性也就越高, 進(jìn)一步進(jìn)行兩兩模塊之間的ME相關(guān)性分析(圖5-C); 將模塊中的所有基因和ME基因在所有樣本中分別進(jìn)行表達(dá)水平分析發(fā)現(xiàn), 各個(gè)模塊內(nèi)的基因其在不同組織部位中的表達(dá)情況呈現(xiàn)不同的特性。每個(gè)模塊內(nèi)的基因表達(dá)水平高度相關(guān), ME表達(dá)水平與模塊整體表達(dá)水平也高度相關(guān)(圖5-D和附圖1), 說明目標(biāo)模塊的ME可充分代表其模塊中的整體基因。

2.4.3 目標(biāo)基因模塊中的關(guān)鍵基因 為獲得與組織部位生長發(fā)育和功能高度相關(guān)的5個(gè)模塊中的核心基因, 利用Cytoscape軟件對(duì)基因互作調(diào)控網(wǎng)絡(luò)進(jìn)行可視化處理, 篩選出模塊中共表達(dá)權(quán)重>0.1并且連通性排名前10的基因, 并結(jié)合NCBI數(shù)據(jù)庫確定模塊內(nèi)的關(guān)鍵基因(圖6)。Turquoise模塊中確定了6個(gè)關(guān)鍵基因, 分別為、、、、、。Green模塊中篩選出了4個(gè)核心基因, 分別為、、、。Magenta模塊中篩選出了4個(gè)核心基因, 分別為、、、。Pink模塊中篩選出了6個(gè)核心基因, 分別為、、、、、。Yellow模塊中篩選出了8個(gè)核心基因, 分別為、、、、、、、。

2.5 候選基因的篩選

檢測到顯著性SNP位點(diǎn)上下游500 kb的LD置信區(qū)間內(nèi)尋找甘藍(lán)型油菜基因, 結(jié)合轉(zhuǎn)錄組差異表達(dá)基因, 通過WGCNA篩選到的與組織部位生長發(fā)育高度相關(guān)的關(guān)鍵模塊中的核心基因(hub genes), 初步確定與甘藍(lán)型油菜生物產(chǎn)量相關(guān)的基因, 將這些基因的蛋白質(zhì)序列與擬南芥基因蛋白質(zhì)序列進(jìn)行BLAST對(duì)比, 結(jié)合前人已報(bào)道的擬南芥同源基因的功能以及特性, 篩選到生物產(chǎn)量性狀中發(fā)揮重要作用的候選基因, 在擬南芥中的同源基因分別是、、、、、(表3)。

和在擬南芥中的同源基因是, 編碼擬南芥葉綠體-1,6-二磷酸酶[40]。豌豆中葉綠體-1,6-二磷酸果糖()干擾載體的轉(zhuǎn)基因植株表現(xiàn)出葉鮮重明顯增加、光合作用的測量結(jié)果顯示具有更高的碳同化率, 表明FBPase作為二氧化碳同化中的關(guān)鍵酶的作用, 并且還可以協(xié)調(diào)碳和氮的代謝[41]。在擬南芥中的同源基因是, 編碼S-腺苷同型半胱氨酸水解酶(SAHH)。SAHH是維持真核生物甲基化穩(wěn)態(tài)的關(guān)鍵酶[42], 會(huì)競爭性地抑制甲基轉(zhuǎn)移酶(MT)活性, 使植株生長緩慢, 繁殖力低, 種子發(fā)芽減少[43-44]。在擬南芥的同源基因是, 編碼在葉綠體基質(zhì)中發(fā)現(xiàn)的一種小肽, 在黑暗中氧化的CP12與甘油醛-3-磷酸脫氫酶和磷酸布洛激酶(碳同化循環(huán)的兩種酶)形成一種非活性超分子復(fù)合物[45],的轉(zhuǎn)錄物在決定成熟葉片和光合作用能力方面起著重要作用[46-48], Marri等人發(fā)現(xiàn), CP12蛋白氧化對(duì)植物不同環(huán)境條件下光合過程的整體平衡至關(guān)重要[49]。在擬南芥中的同源基因是一種卡爾文循環(huán)酶, 具有光合固碳作用[48]。在擬南芥中的同源基因?yàn)? 磷酸化甘油醛-3-P脫氫酶(GAPC-1)是一種高度保守的胞漿酶, 純合缺失突變株表現(xiàn)出生長延遲、角果形態(tài)改變和種子數(shù)量低[51]。在擬南芥中的同源基因是, 是擬南芥營養(yǎng)組織中主要肌動(dòng)蛋白基因, 對(duì)植物生長發(fā)育有著極其重要的作用[52]。多項(xiàng)研究結(jié)果表明, 甘藍(lán)型油菜生物產(chǎn)量相關(guān)候選基因的同源基因功能, 在植物氧化還原、能量和碳水化合物代謝、光合作用、物質(zhì)積累等生長發(fā)育過程中發(fā)揮作用。

A: 樣本聚類與數(shù)據(jù)矯正; B和C: 基于規(guī)模獨(dú)立性和均值連通性選擇軟閾值; D: 已識(shí)別模塊的樹狀聚類圖; E: 已識(shí)別模塊的熱圖。

A: sample clustering to detect outliers; B and C: the selection of the soft threshold based on scale independence and mean connectivity; D: cluster dendrogram of the identified modules; E: the heatmap of identified modules.

A: 基因共表達(dá)網(wǎng)絡(luò)模塊與組織部位關(guān)聯(lián)熱圖; B: 各樣本中Turquoise模塊的所有基因與相應(yīng)ME的表達(dá)水平; C: 不同模塊兩兩之間ME的相關(guān)性; D: ME聚類樹。

A: association analysis of gene co-expression network modules with tissues; B: expression levels of all genes and corresponding ME in turquoise module of each sample; C: ME correlation between different modules; D: ME cluster tree.

(圖6)

(圖6)

表3 生物產(chǎn)量候選基因

3 討論

生物產(chǎn)量體現(xiàn)作物總的生產(chǎn)力, 對(duì)作物的經(jīng)濟(jì)產(chǎn)量有著密不可分的關(guān)系。本研究所使用甘藍(lán)型油菜自然群體中的生物產(chǎn)量變化范圍廣、變異系數(shù)大, 說明生物產(chǎn)量是比較復(fù)雜的性狀。本研究在2017、2019年和BLUP環(huán)境下共檢測到9個(gè)生物產(chǎn)量顯著關(guān)聯(lián)標(biāo)記位點(diǎn), 貢獻(xiàn)率在5.64%~7.98%。與生物產(chǎn)量顯著關(guān)聯(lián)的SNP位點(diǎn)分布在各個(gè)染色體, Lu等[53]利用520份來自世界各地的甘藍(lán)型油菜材料作為自然群體, 使用MLM進(jìn)行全基因組關(guān)聯(lián)分析, 在除了A01、A06、C02、C07以外的15條染色體上共檢測到26個(gè)生物產(chǎn)量SNP, 貢獻(xiàn)率在4.33%~17.03%之間; 60K SNP芯片對(duì)520分材料構(gòu)成的自然群體進(jìn)行全基因組關(guān)聯(lián)分析, 在除了A01、A02和A10意外的16條染色體共檢測到89個(gè)生物產(chǎn)量SNP位點(diǎn), 貢獻(xiàn)率在3.47%~8.36%之間[54]; Luo等[55]利用155分甘藍(lán)型油菜材料構(gòu)成的自然群體并使用60K SNP芯片檢測到1個(gè)生物產(chǎn)量SNP, 其貢獻(xiàn)率為0.76%。本研究中, 生物產(chǎn)量在2017、2019年和BLUP環(huán)境下共檢測到9個(gè)顯著關(guān)聯(lián)標(biāo)記位點(diǎn), 分別位于A03、A07、C01、C03、C04、C06、C09染色體, 貢獻(xiàn)率在5.64%~7.98%。前人檢測到與生物產(chǎn)量顯著關(guān)聯(lián)的部分SNP位點(diǎn)也定位在這些染色體上, 說明生物產(chǎn)量的構(gòu)成因素比較復(fù)雜, 且由多個(gè)SNP位點(diǎn)協(xié)同控制。

通過RNA-Seq對(duì)極端表型材料的差異表達(dá)基因分析表明, 葉片和角果皮主要富集到光合作用、淀粉和蔗糖代謝等生物途徑, 這些代謝途徑可能在油菜植株的物質(zhì)積累過程中發(fā)揮重要作用, 因此, 這些差異基因的不同表達(dá)模式可能是造成2個(gè)材料生物產(chǎn)量有顯著差異的主要原因。莖稈中的差異表達(dá)基因主要富集在細(xì)胞壁合成等途徑, 油菜莖稈次生細(xì)胞壁的發(fā)育對(duì)其抗倒伏有著決定性影響, 莖稈的發(fā)育狀況良好對(duì)油菜碳水化合物的運(yùn)輸積累有著積極作用[55]。因此, 本研究極端材料莖稈中的差異基因有可能在這些方面影響其生物產(chǎn)量。種子是油菜主要的儲(chǔ)存和收獲器官, 其發(fā)育狀況直接影響著最終產(chǎn)量[56], 因此也直接影響生物產(chǎn)量。本研究極端生物產(chǎn)量材料種子中的差異基因主要富集在營養(yǎng)庫活性、蔗糖淀粉代謝等功能和途徑。

WGCNA是共表達(dá)網(wǎng)絡(luò)分析有效的分析方法, 能夠特異地篩選出與目標(biāo)性狀具有高度生物學(xué)意義的共表達(dá)模塊, 在玉米等植物中已經(jīng)被證明是一種高效的數(shù)據(jù)挖掘方法。本研究得到15個(gè)生物產(chǎn)量加權(quán)共表達(dá)網(wǎng)絡(luò)模塊, 分析發(fā)現(xiàn)一些模塊與部分組織有著顯著的相關(guān)性, 分別是葉片與Turquoise模塊, 莖稈與Green、Magenta、Pink模塊, 開花后30 d主軸種子與Yellow模塊相關(guān)。葉片、莖稈和種子均在作物的生物產(chǎn)量的構(gòu)成中發(fā)揮極其重要的作用, 直接影響到農(nóng)產(chǎn)品最終的產(chǎn)量。結(jié)合GWAS和轉(zhuǎn)錄組差異基因得到的結(jié)果、整合WGCNA挖掘得到關(guān)聯(lián)模塊中的核心基因(hub genes), 最終篩選到一批可能與甘藍(lán)型油菜生物產(chǎn)量密切相關(guān)的基因, 分別是、、、、、、。這些基因的干擾、敲除或者超量表達(dá)均引起了植株的整體重量的變化[40-52], 在擬南芥生物產(chǎn)量中發(fā)揮重要作用。

由于生物組學(xué)大數(shù)據(jù)存在復(fù)雜、多層次和信息互補(bǔ)的特點(diǎn), 分析這些數(shù)據(jù)的一個(gè)關(guān)鍵目標(biāo)是確定可預(yù)測表型性狀的有效模型, 發(fā)現(xiàn)重要性狀的關(guān)鍵調(diào)節(jié)因子并闡明其生物功能。本研究根據(jù)GWAS分析得到的與生物產(chǎn)量高度關(guān)聯(lián)SNP標(biāo)記, 同時(shí)結(jié)合轉(zhuǎn)錄組分析確定的差異基因以及WGCNA得到的基因富集模塊進(jìn)一步篩選出與甘藍(lán)型油菜生物產(chǎn)量相關(guān)的基因。通過BLAST對(duì)比找到這些基因在擬南芥中的同源基因, 并根據(jù)其注釋以及已有報(bào)道進(jìn)一步篩選與生物產(chǎn)量高度相關(guān)的候選基因, 為油菜中生物產(chǎn)量候選基因的功能研究奠定基礎(chǔ)。

4 結(jié)論

本研究通過2年間甘藍(lán)型油菜自然群體中生物產(chǎn)量對(duì)大多數(shù)產(chǎn)量相關(guān)性狀都具有正向效應(yīng), 說明甘藍(lán)型油菜的生物產(chǎn)量是其他產(chǎn)量相關(guān)性狀的基礎(chǔ)和保障。結(jié)合GWAS關(guān)聯(lián)基因與轉(zhuǎn)錄組差異基因得到178個(gè)基因, 根據(jù)擬南芥同源基因的功能, 篩選得到與生物產(chǎn)量性狀相關(guān)的和為關(guān)鍵候選基因, 它們?cè)谀芰亢吞妓衔锎x以及光合作用等方面中發(fā)揮作用。此外, 選取WGCNA所得到的重點(diǎn)模塊中連通性前10的基因作為hub gene, 通過同源基因功能注釋分析, 得到生物產(chǎn)量候選基因、、、和, 在光合作用的卡爾文循環(huán)、碳同化、物質(zhì)積累等方面發(fā)揮作用。

[1] 金森森, 羅帥. 全球視角下的南美大豆貿(mào)易發(fā)展趨勢(shì). 中國糧食經(jīng)濟(jì), 2014, (4): 30–32.Jin S S, Luo S. The development trend of South American soybean trade from a global perspective., 2014, (4): 30–32 (in Chinese with English abstract).

[2] 陶伯玉. 甘藍(lán)型油菜主要農(nóng)藝、品質(zhì)性狀的遺傳分析和QTL定位. 南京農(nóng)業(yè)大學(xué)碩士學(xué)位論文, 江蘇南京, 2015. Tao B Y. Genetic Analysis and QTL Mapping of Main Agronomic and Quality Traits in. MS Thesis of Nanjing Agricultural University, Nanjing, Jiangsu, China, 2015 (in Chinese with English abstract).

[3] 張衛(wèi). 明確大宗油料作物發(fā)展規(guī)劃提升國內(nèi)大宗油料作物自給力. 中國食品, 2016, (19): 82–89. Zhang W. Clarify the development plan of bulk oil crops to enhance the self-sufficiency of domestic bulk oil crops., 2016, (19): 82–89 (in Chinese with English abstract).

[4] 徐露萍. 淺談?dòng)筒朔N植技術(shù)的應(yīng)用與推廣. 南方農(nóng)業(yè), 2015, 9(3): 17–18. Xu L P. Talking about the application and extension of rape planting technology., 2015, 9(3): 17–18 (in Chinese with English abstract).

[5] 李培武, 丁小霞, 張文, 謝立華, 周海燕. 中加雙低油菜質(zhì)量標(biāo)準(zhǔn)體系比較及我國發(fā)展對(duì)策探討. 農(nóng)產(chǎn)品質(zhì)量與安全, 2006, (6): 23–26. Li P W, Ding X X, Zhang W, Xie L H, Zhou H Y. Comparison of quality standard systems for double-low rapeseed between China and Canada and discussion on my country’s development countermeasures., 2006, (6): 23–26 (in Chinese with English abstract).

[6] Mirza H, Karim F M. Rapeseed () cultivation in dry season., 2010, 12: 27–80.

[7] 傅廷棟. 油菜的品種改良. 作物研究, 2007, 21(3): 159–162. Fu T D. Rapeseed variety improvement., 2007, 21(3): 159–162 (in Chinese with English abstract).

[8] Liersch A, Bocianowski J, Henryk W, Laurencja S. Assessment of genetic relationships in breeding lines and cultivars ofand their implications for breeding winter oilseed rape., 2016, 56: 15–40.

[9] Zhao J, Meng J. Genetic analysis of loci associated with partial resistance toin rapeseed (L.)., 2003, 106: 759–764.

[10] Wen Y C, Fu T D, Tu J X, MA C Z, Shen J X, Wen J, Zhang S F. Factors analysis of silique shatter resistance in rapeseed (L.)., 2010, 32: 25–29.

[11] Chen S, Qian T. Analysis on technical efficiency and influencing factors of rapeseed production., 2016, 8: 1–7.

[12] Jiang H D, Zhou Q, Li N, Sun X F. Effect of Cd on the growth and physiological characteristics of rape seedlings., 2006, 28: 39–43.

[13] 張哲, 殷艷, 劉芳, 王積軍, 傅廷棟. 我國油菜多功能開發(fā)利用現(xiàn)狀及發(fā)展對(duì)策. 中國油料作物學(xué)報(bào), 2018, 40: 618–623. Zhang Z, Yin Y, Liu F, Wang J J, Fu T D. The current status and development countermeasures of multifunctional development and utilization of rapeseed in my country., 2018, 40: 618–623 (in Chinese with English abstract).

[14] Chen H, He H, Zou Y J, Chen W, Yu R B, Liu X, Yang Y, Gao Y M, Xu J L, Fan L M, Li Y, Li Z K, Deng X W. Development and application of a set of breeder-friendly SNP markers for genetic analyses and molecular breeding of rice (L.)., 2011, 123: 869–879.

[15] Ma T Q, Lin L B. Application of molecular marker technologies in rape breeding., 2004, 2: 728–732.

[16] Xiong Q F, Wen J, Li X H, Shen J X. Technological innovation and industrial development of rapeseed in China., 2014, 16: 14–22.

[17] Beecher H G, Dunn B W, Thompson J A, Humphreys E, Mathews S K, Timsina J. Effect of raised beds, irrigation and nitrogen management on growth, water use and yield of rice in south- eastern Australia., 2006, 46: 1363–1372.

[18] Xu Z, Yu Z, Zhao J. Theory and application for the promotion of wheat production in China: past, present and future., 2013, 93: 2339–2350.

[19] Cheng L, Li H P, Qu B, Huang T, Tu J X, Fu T D, Liao Y C. Chloroplast transformation of rapeseed () by particle bombardment of cotyledons., 2010, 29: 371–381.

[20] Wang G C, Liu Z, Yang G S. Development of high oil content lines ofthrough isolated microspore culture., 2007, 29: 382–386.

[21] Liao X, Wang H Z. Present status of rapeseed science and technology innovation system in China and suggestion for further development., 2007, 6: 28–34.

[22] Cheng L, Li H P, Qu B. Chloroplast transformation of rapeseed () by particle bombardment of cotyledons., 2010, 29: 371–381.

[23] Zhang X K, Chen J, Chen L, Wang H Z, Li J N. Imbibition behavior and flooding tolerance of rapeseed seed (L.) with different testa color., 2008, 55: 1175–1184.

[24] Palmieri S, Leoni O. Industrial prospectives of rapeseed oil. Physiologic technological and qualitative aspects., 1992, 14: 63–70.

[25] Weilenmann M E, Lúquez J, Suárez W. Variability of biological yield, economic yield and harvest index in soybean cultivars grown in irrigated and rainfed environments in Argentina., 2000, 21: 39–40.

[26] Li S Y, Zhu R K, Varshney J, Zhan X, Zheng J, Shi X, Wang G, Wang H. A systematic dissection of the mechanisms underlying the natural variation of silique number in rapeseed (L.) germplasm., 2020, 18: 568–580.

[27] Luo X, Ma C Z, Yue Y, Hu K N, Li Y Y, Duan Z Q, Wu M, Tu J X, Shen J X, Yi B, Fu T D. Unravelling the complex trait of harvest index in rapeseed (L.) with association mapping., 2015, 16: 379–388.

[28] 姚東和, 楊民勝, 李志輝. 林分密度對(duì)巨尾桉生物產(chǎn)量及生產(chǎn)力的影響. 中南林業(yè)科技大學(xué)學(xué)報(bào), 2000, 20(3): 20–23. Yao D H, Yang M S, Li Z H. The effect of stand density on the biological yield and productivity of., 2000, 20(3): 20–23(in Chinese with English abstract).

[29] 李躍建, 宋荷仙, 朱華忠. 小麥?zhǔn)斋@指數(shù)、生物產(chǎn)量和籽粒產(chǎn)量的穩(wěn)定性分析. 西南農(nóng)業(yè)學(xué)報(bào), 1998, (1): 25–30. Li Y J, Song H X, Zhu H Z. Analysis of stability of wheat harvest index, biological yield and grain yield., 1998, (1): 25–30 (in Chinese with English abstract).

[30] Dai X L, Zhao J X, Xiang Y, Zhang T, Ren T B, Cheng G P. Correlation analysis of harvest index and plant traits of hybrid., 2018, 2: 47–58.

[31] Tian J, Ma K, Saaem I. Advancing high-throughput gene synthesis technology., 2009, 5: 714–722.

[32] Dillies M, Rau A, Aubert J, Hennequet C, Jeanmougin M, Servant A, Keime C, Marot G, Castel D, Estelle J, Guernec G, Jagla B, Jouneau L, Lalo? C, Gall L, Scha?ffer B, Le C, Guedj M, Jaffrézic F, Consortium. A comprehensive evaluation of normalization methods for Illumina high-throughput RNA sequencing data analysis., 2013, 14: 671–683.

[33] Wang J X, Zhang X Q, Shi M L, Gao L J, Niu X F, Te R, Chen L, Zhang W W. Metabolomic analysis of the salt-sensitive mutants reveals changes in amino acid and fatty acid composition important to long-term salt stress insp. PCC 6803., 2014, 14: 431–440.

[34] 秦天元, 孫超, 畢真真, 梁文君, 李鵬程, 張俊蓮, 白江平. 基于WGCNA的馬鈴薯根系抗旱相關(guān)共表達(dá)模塊鑒定和核心基因發(fā)掘. 作物學(xué)報(bào), 2020, 46: 1033–1051. Qin T Y, Sun C, Bi Z Z, Liang W J, Li P C, Zhang J L, Bai J P. Identification of drought-related co-expression modules and hub genes in potato roots based on WGCNA., 2020, 46: 1033–1051 (in Chinese with English abstract).

[35] 宋長新, 雷萍, 王婷. 基于WGCNA算法的基因共表達(dá)網(wǎng)絡(luò)構(gòu)建理論及其R軟件實(shí)現(xiàn). 基因組學(xué)與應(yīng)用生物學(xué), 2013, 32: 135–141. Song C X, Lei P, Wang T. Gene co-expression network construction theory based on WGCNA algorithm and its R software implementation., 2013, 32: 135–141 (in Chinese with English abstract).

[36] Lu K, Wei L J, Li X L, Wang Y T, Wu J, Liu M, Zhang C, Chen Z Y, Xiao Z C, Jian H J, Cheng F, Zhang K, Du H, Cheng X C, Qu C M, Qian W, Liu L Z, Wang R, Zou Q Y, Ying J M, Xu X F, Mei J Q, Liang Y, Chai Y R, Tang Z L, Wan H F, Ni Y, He Y J, Lin N, Fan Y H, Sun W, Li N N, Zhou G, Zheng H K, Wang X W, Andrew H P, Li J N. Whole-genome resequencing revealsorigin and genetic loci involved in its improvement., 2019, 10: 1154–1165.

[37] Bradbury P J, Zhang Z, Kroon D E, Casstevens T M, Ramdoss Y, Buckler E S. TASSEL: software for association mapping of complex traits in diverse samples., 2007, 23: 2633–2635.

[38] Barrett J C, Fry B, Maller J, Daly M J. Haploview: analysis and visualization of LD and haplotype maps., 2005, 21: 263–265.

[39] Kohl M, Wiese S, Warscheid B. Cytoscape: software for visualization and analysis of biological networks., 2011, 696: 291–303.

[40] Livingston A K, Cruz J A, Kohzuma K, Dhingra A. Anmutant with high cyclic electron flow around photosystem I (hcef) involving the NADPH dehydrogenase complex., 2010, 22: 221–233.

[41] Mariam S. Increased sucrose level and altered nitrogen metabolism intransgenic plants expressing antisense chloroplastic fructose-1,6-bisphosphatase.,2004, 55: 2495–2503.

[42] Mull L, Ebbs M L, Bender J. A histone methylation-dependent DNA methylation pathway is uniquely impaired by deficiency inS-adenosylhomocysteine hydrolase., 2005, 174: 1161–1171.

[43] Ouyang B, Joung J G, Kolenovsky A, Koh C, Nowak J, Caplan A, Keller W A, Cui Y, Cutler A J, Tsang E W. Transcriptome profiling and methyl homeostasis of anmutant deficient in S-adenosylhomocysteine hydrolase1 (SAHH1)., 2012, 79: 315–331.

[44] Godge M R, Kumar D, Kumar P P.HOG1 gene and its petunia homolog PETCBP act as key regulators of yield parameters., 2008, 27: 1497–1507.

[45] Pedro S C F R, Mazhar S, Rosalba M, Mathilde F, Stéphanie B, Rebecca L, Barbara M, Conrad W, Hervé V, Ian F. TheH OMOLOGY-DEPENDENT GENE SILENCING1 gene codes for an-adenosyl--homocysteine hydrolase required for DNA methylation-dependent gene silencing., 2005, 17: 404–417.

[46] López C, Abuzaid O, Lawson T, Raines C A.CP12 mutants have reduced levels of phosphoribulokinase and impaired function of the calvin-benson cycle., 2017, 68: 2285–2298.

[47] Fermani S, Trivelli X, Sparla F, Thumiger A, Calvaresi M, Marri L, Falini G, Zerbetto F, Trost P. Conformational selection and folding-upon-binding of intrinsically disordered protein CP12 regulate photosynthetic enzymes assembly., 2012, 287: 21372–21383.

[48] Singh P, Kaloudas D, Raines C A. Expression analysis of theCP12 gene family suggests novel roles for these proteins in roots and floral tissues.,2008, 59: 3975–3985.

[49] Marri L, Thieulin P G, Lebrun R, Puppo R, Zaffagnini M, Trost P, Gontero B, Sparla F. CP12-mediated protection of Calvin-Benson cycle enzymes from oxidative stress., 2014, 97: 228–237.

[50] Ferro M, Salvi D, Seigneurin B D, Court M, Moyet L, Ramus C, Miras S, Mellal M, Le G S, KiefferJ S, Bruley C, Garin J, Joyard J, Masselon C, Rolland N. AT_CHLORO, a comprehensive chloroplast proteome database with subplastidial localization and curated information on envelope proteins.,2010, 9: 1063–1084.

[51] Henry E, Fung J, Liu G, Drakakak G, Coaker G. Beyond glycolysis: GAPDHs are multi-functional enzymes involved in regulation of ROS, autophagy, and plant immune responses., 2015, 11: e1005199.

[52] Ringli C, Baumberger N, Diet A, Frey B, Keller B. ACTIN2 is essential for bulge site selection and tip growth during root hair development of., 2002, 129: 1464–1472.

[53] Lu K, Liu P, Zhang C, Lu J H, Yang B, Xiao Z C, Liang Y, Xu X F, Qu C M, Zhang K, Liu L Z, Zhu Q L, Fu M L, Yuan X Y, Li J N. Genome-wide association and transcriptome analyses reveal candidate genes underlying yield-determining traits in., 2017, 8: 206–220.

[54] Lu K, Xiao Z C, Jian H J, Liu P, Qu C M, Fu M L, He B, Tie L M, Liang Y, Xu X F, Li J N. A combination of genome-wide association and transcriptome analysis reveals candidate genes contro-lling harvest index-related traits in., 2016, 6: 36452.

[55] Sun Y Y, Liu T T, Yang H Y, Zuo Q S, Zhou G S, Wu J S. The correlation between the growth characteristics of rapeseed stalk and its yield and the lodging resistance., 2014, 53: 30–35.

[56] 柳寒. 油菜灌漿期種子和角果皮基因表達(dá)差異及相關(guān)基因的功能分析. 華中農(nóng)業(yè)大學(xué)博士學(xué)位論文, 湖北武漢, 2015.Liu H. Gene Expression and Functional Analysis in Seeds and Siliques of Rapeseed during the Filling Stage. PhD Dissertation of Huazhong Agricultural University, Wuhan, Hubei, China, 2015 (in Chinese with English abstract).

附圖1 各樣本中關(guān)鍵模塊的所有基因與相應(yīng)ME的表達(dá)水平

Fig. S1 Expression levels of all genes and corresponding ME in key modules of each sample

Integrating GWAS and WGCNA to screen and identify candidate genes for biological yield inL.

WANG Yan-Hua1,2,**, LIU Jing-Sen1,2,**, and LI Jia-Na1,2,*

1College of Agronomy and Biotechnology, Southwest University / Chongqing Engineering Research Center for Rapeseed, Chongqing 400715, China;2Academy of Agricultural Sciences, Southwest University, Chongqing 400715, China

Biomass yield is especially important for, as it is the basis for high yields of crops. In this study, the phenotypic data of the natural populations composed of 588 materials were used for genome-wide association analysis (GWAS). We performed the transcriptome sequencing (RNA-seq) of biomass yield using ‘CQ45’ (high biological yield material) and ‘CQ46’ (low biological yield material). A weighted gene co-expression network analysis (WGCNA) network was constructed by integrating transcriptome data of six tissues of the extreme materials, such as stalks, leaves, 30 day after flowering (DAF) seeds of main inflorescence and lateral branch, 30 DAF pod keratin of main branch and lateral branch. We finally screened the candidate genes related to biomass yield. The main results are as follows: Biomass yields inhad positive effects on most yield-related traits; K + PCA model was the best model for biomass analysis of the natural population, and nine significant loci were detected in the best model (< 1/385691 or< 0.05/385691); according to 36 groups of transcriptome data, MAD value of each gene was calculated. A total of 5052 genes with MAD value of the top 5% were selected to construct WGCNA. Fifteen gene modules were obtained, among which, five genes co-expression modules were significantly correlated with leaves, stems, and seeds of 30 DAF. The hub genes of the key modules in WGCNA, the significant SNP loci obtained from GWAS, and the extreme phenotypic differential genes were integrated to identify the candidate genes. Theirhomologous genes were,,, and, which played the important roles in the Calvin cycle, carbon assimilation, and material accumulation of photosynthesis.

GWAS; WGCNA; RNA seq;

10.3724/SP.J.1006.2021.04175

本研究由中國博士后科學(xué)基金面上項(xiàng)目(2019M653319), 重慶市自然科學(xué)基金博士后科學(xué)基金項(xiàng)目(cstc2019jcyj-bshXO116)和高等學(xué)校學(xué)科創(chuàng)新引智基地項(xiàng)目(“111”項(xiàng)目)(B12006)資助。

This study was supported by the Project of China Postdoctoral Science Foundation (2019M653319), the Project of Chongqing Natural Science Foundation Postdoctoral Science Foundation (cstc2019jcyj-bshXO116), and the Project of Intellectual Base for Discipline Innovation in Colleges and Universities (“111” Project) (B12006).

李加納, E-mail: ljn1950@swu.edu.cn

**同等貢獻(xiàn)(Contributed equally to this work)

E-mail: hawer313@163.com

2020-07-30;

2020-12-01;

2021-01-11.

URL: https://kns.cnki.net/kcms/detail/11.1809.S.20210108.1700.010.html