大豆紅色種皮的色素鑒定和基因定位

曹杰，谷勇哲，洪慧龍，吳海濤，張霞，孫建強，包立高，邱麗娟

大豆紅色種皮的色素鑒定和基因定位

曹杰1,2，谷勇哲2，洪慧龍2，吳海濤2，張霞2，孫建強3，包立高4，邱麗娟1,2

1吉林農(nóng)業(yè)大學(xué)生命科學(xué)學(xué)院，長春 130118；2中國農(nóng)業(yè)科學(xué)院作物科學(xué)研究所，北京 100081；3東北農(nóng)業(yè)大學(xué)農(nóng)學(xué)院，哈爾濱 150030；4內(nèi)蒙古自治區(qū)農(nóng)牧業(yè)技術(shù)推廣中心，呼和浩特 010018

【目的】揭示種子發(fā)育過程中種皮花青素（anthocyanin）的含量變化以及導(dǎo)致泰興矮腳紅（TXAJH）紅色種皮的主要花青素成分；定位控制花青素合成積累的關(guān)鍵基因，為深入了解紅色種皮形成的調(diào)控機制奠定基礎(chǔ)?！痉椒ā坷贸咝б合嗌V串聯(lián)質(zhì)譜（ultra-high performance liquid chromatography-tandem mass spectrometry，UPLC-ESI-MS/MS）檢測黃色種皮大豆綏農(nóng)14（SN14）和紅色種皮大豆TXAJH不同發(fā)育階段種皮的花青素成分與含量，分析與種皮顏色變化密切相關(guān)的花青素成分；利用SN14和TXAJH雜交構(gòu)建的重組自交系（recombinant inbred lines，RIL）群體進(jìn)行分離群體分組混合分析（bulked segregant analysis，BSA），初步定位紅色種皮相關(guān)基因的候選區(qū)域，在此基礎(chǔ)上，結(jié)合標(biāo)記連鎖分析縮小候選區(qū)間并預(yù)測紅色種皮候選基因；最后通過qRT-PCR驗證候選基因的表達(dá)情況?！窘Y(jié)果】檢測SN14和TXAJH 4個發(fā)育階段的種皮，共發(fā)現(xiàn)12種花青素。在成分上，總花青素的聚類分析表明，TXAJH與SN14之間以及TXAJH顯色前后之間的種皮花青素組成均存在明顯差異。在含量上，種子發(fā)育過程中，SN14種皮花青素的含量逐漸下降，而TXAJH種皮的含量迅速升高并保持穩(wěn)定，種皮顯色后，二者的花青素含量呈現(xiàn)極顯著差異，在成熟階段，TXAJH種皮花青素的含量是SN14的200倍以上。矢車菊素-3-O-葡萄糖苷（Cyanidin-3-O-glucoside，Cy-3-glu）、芍藥花素-3-O-葡萄糖苷（Peonidin-3-O-glucoside，Pn-3-glu）和牽?；ㄋ?3-O-葡萄糖苷（Petunidin-3-O-glucoside，Pt-3-glu）是導(dǎo)致TXAJH種皮呈現(xiàn)紅色的重要原因。BSA-seq關(guān)聯(lián)分析將紅色種皮基因的候選區(qū)間定位于第8染色體上，長度為8.66 Mb。利用27個多態(tài)性標(biāo)記進(jìn)行連鎖分析得到10種單倍型，最終將候選區(qū)間縮小至702 kb。該區(qū)間中在親本間存在非同義變異的基因共37個，其中，編碼MYB轉(zhuǎn)錄因子，和編碼bHLH轉(zhuǎn)錄因子，它們可能參與花青素的生物合成調(diào)控；編碼花青素還原酶1，可以將花青素轉(zhuǎn)化為原花青素（proanthocyanidin，PA）?；虮磉_(dá)分析結(jié)果表明，候選基因和花青素生物合成途徑相關(guān)基因在SN14與TXAJH中的表達(dá)模式相似，均為前者低于后者。種皮花青素主要成分與候選基因表達(dá)水平的關(guān)聯(lián)分析結(jié)果顯示二者之間存在極強的相關(guān)性?！窘Y(jié)論】SN14與TXAJH的種皮花青素組成存在差異，TXAJH紅色種皮呈現(xiàn)紅色可能是Cy-3-glu、Pn-3-glu和Pt-3-glu積累的結(jié)果。預(yù)測、、和為紅色種皮候選基因，其中、和可能對花青素生物合成途徑的多個基因產(chǎn)生調(diào)控作用。

大豆；種皮色；花青素；BSA-seq；基因定位；轉(zhuǎn)錄因子

0 引言

【研究意義】大豆（(L.) Merr.）是重要的經(jīng)濟作物，可以為人類和牲畜提供優(yōu)質(zhì)的植物蛋白和油脂[1-2]。大豆種皮顏色是一類非常便于觀察的生物學(xué)性狀，在遺傳學(xué)研究中常作為形態(tài)標(biāo)記[3]。使用不同顏色來標(biāo)記區(qū)分用于制藥、食品飼料加工等商業(yè)用途的大豆將有望極大地降低原料篩選成本[4]。植物種皮的顏色主要取決于所合成花青素（anthocyanin）的種類和含量[5]。花青素是具有良好抗氧化活性的天然色素，長期食用花青素或富含花青素的食物對健康有一系列的好處，包括但不限于保護(hù)心血管、保護(hù)神經(jīng)、改善視力、促進(jìn)新陳代謝、抗菌抗炎癥以及預(yù)防和治療癌癥等[6-12]。此外，花青素可以有效清除植物中由于非生物脅迫而產(chǎn)生的過量活性氧（reactive oxygen species，ROS），降低干旱、低溫、高鹽和重金屬污染等逆境對植物造成的直接或間接損傷[13-17]。大豆中紅色種皮較為罕見，鑒定其中的花青素成分，發(fā)掘與紅色種皮形成相關(guān)基因，對解析種皮顏色形成機制，提升大豆品質(zhì)和經(jīng)濟效益具有一定的參考意義?！厩叭搜芯窟M(jìn)展】經(jīng)典遺傳研究表明，大豆籽粒顏色是一個受多位點控制的復(fù)雜性狀，已發(fā)現(xiàn)至少9個位點（、、、、、、、和）參與控制大豆種皮顏色[18]。、和位點通過影響葉綠素代謝，使種皮表現(xiàn)綠色[3, 19]。、、、、和位點通過控制大豆種皮中花青素的合成與分布使種皮的整體或局部發(fā)生色素沉積進(jìn)而形成各種顏色的種皮[20-28]?；ㄇ嗨厥穷慄S酮（flavonoids）物質(zhì)下面的重要子類，在可見光波段具有廣泛的吸收范圍，是使植物產(chǎn)生顏色的重要物質(zhì)，通常導(dǎo)致種皮呈現(xiàn)棕色、紅色、雙色、黑色等顏色[5]。在影響大豆種皮顏色的遺傳位點中，是查爾酮合成酶（chalcone synthase，CHS）基因的重復(fù)區(qū)域[20]。CHS是催化花青素生物合成途徑初始步驟的專用酶，多個的存在會降低其自身在種皮中的表達(dá)水平進(jìn)而抑制種皮的花青素合成[21]。和位點分別與黃烷酮3′-羥化酶（flavonoid 3′-hydroxylase，F(xiàn)3′H）基因和黃烷酮3′5′-羥化酶（flavonoid 3′5′-hydroxylase，F(xiàn)3′5′H）基因共分離[22-24]。F3′H和F3′5′H對黃烷酮B環(huán)不同位置的羥基化修飾將決定最終合成花青素的顏色類型，對植物顏色的呈現(xiàn)有重要影響[25-26]。和位點分別與編碼R2R3-MYB轉(zhuǎn)錄因子和Argonaute5（AGO5）蛋白的基因有關(guān)[27-28]。R2R3-MYB轉(zhuǎn)錄因子調(diào)控花青素合成途徑部分關(guān)鍵基因的表達(dá)，可導(dǎo)致純黑/純褐色種皮或雙色條紋種皮[27]。AGO5蛋白具有調(diào)控siRNAs在種皮中分布的功能，使種皮特定部位例如種臍和鞍區(qū)發(fā)生色素沉積形成鞍掛種皮[28]。位點與紅棕色種皮有關(guān)，但只在和位點均為隱性的時候才影響種皮顏色[18]。在基因定位方法上，傳統(tǒng)的QTL圖譜通常涉及使用分布于整個基因組中的分子標(biāo)記對群體中的大量個體進(jìn)行基因型和表型分析以保證足夠的統(tǒng)計能力，這一過程需要耗費大量的時間和人力物力成本[29]。相比之下，利用分離體分組混合分析（bulked segregant analysis，BSA）方法可以通過構(gòu)建目標(biāo)性狀極端差異的基因混池快速篩選與目標(biāo)性狀緊密相關(guān)的位點來定位基因[30-31]。新一代測序技術(shù)的發(fā)展迭代使更多物種陸續(xù)完成全基因組測序，基于全基因組測序的BSA方法可以在沒有遺傳圖譜的情況下，對目標(biāo)性狀進(jìn)行快速精準(zhǔn)的定位分析，具有經(jīng)濟、方便、快捷的優(yōu)點[32-33]。目前，BSA-seq方法已被廣泛應(yīng)用于小麥[34]、水稻[35]、玉米[36]、花生[33, 37]、大豆等的基因定位中，并取得顯著成果?！颈狙芯壳腥朦c】目前，大豆種皮花青素的研究多集中于黑色種皮的成熟種子，而對紅色種皮以及種皮發(fā)育過程花青素含量變化的研究較為罕見，且已發(fā)現(xiàn)的位點/基因尚不能完全解釋種皮顏色的形成機制?！緮M解決的關(guān)鍵問題】本研究通過對紅、黃2種顏色大豆不同發(fā)育階段的種皮進(jìn)行花青素成分含量檢測，揭示種子發(fā)育過程中種皮花青素的含量變化規(guī)律以及決定紅色種皮的主要花青素成分；利用綏農(nóng)14與泰興矮腳紅雜交構(gòu)建的重組自交系（recombinant inbred lines，RIL）群體進(jìn)行BSA-seq和基因定位，預(yù)測紅色種皮顏色候選基因，為深入解析種皮顏色形成的調(diào)控機制奠定基礎(chǔ)。

1 材料與方法

1.1 植物材料

前期以黃色種皮栽培品種綏農(nóng)14（SN14，ZDD22648）為母本，紅色種皮地方品種泰興矮腳紅（TXAJH，ZDD04430）為父本，通過雜交和連續(xù)自交的方法構(gòu)建了F9RIL群體。

SN14和TXAJH以及RIL群體的188個品系于2022年6月在中國農(nóng)業(yè)科學(xué)院作物科學(xué)研究所試驗田（116°20′E，39°57′N）種植，采用常規(guī)田間管理模式。播種后第30天收集親本及每個品系的幼嫩葉片用于BSA-seq和基因定位。分別在開花后第30、40、50、60和70天剝?nèi)N14和TXAJH的種皮用于花青素含量檢測和qRT-PCR分析，每個時期的每個樣品均包含3個生物重復(fù)。以上所有樣品離體后立即在液氮中速凍，-80 ℃保存?zhèn)溆谩?/p>

1.2 種皮花青素提取

分別選取SN14與TXAJH開花后第40、50、60和70天的種皮進(jìn)行花青素成分檢測。樣品經(jīng)充分研磨后，每份取0.1 g凍干粉末溶解于1 mL 70%甲醇中，充分渦旋混勻，4 °C過夜萃取，然后12 000 r/min離心10 min，吸取上清，用微孔濾膜（0.22 μm pore size）過濾得花青素提取液，避光4 ℃保存?zhèn)溆谩?/p>

1.3 花青素檢測

用超高效液相色譜串聯(lián)質(zhì)譜（UPLC-ESI-MS/MS）進(jìn)行花青素檢測。將所有待測樣品的等量混合物作為QC樣本，每8個待測樣本中插入一個QC樣本以監(jiān)測在相同的處理方法下分析過程的重復(fù)性。液相條件主要包括：1）色譜柱：Waters Acquity UPLC HSS T3 C18 1.8 μm，2.1 mm×100 mm；2）流動相：各含0.1%甲酸的超純水（A相）和乙腈（B相）；3）洗脫梯度（v/v）：0 min水/乙腈為95/5，10.0 min為5/95，11.0 min為5/95，11.1 min為95/5，15.0 min為95/5；4）流速0.4 mL·min-1；柱溫40 ℃；進(jìn)樣量2 μL。質(zhì)譜條件主要包括：電噴霧離子源（electrospray ionization，ESI）溫度為550 ℃，質(zhì)譜電壓為5 500 V/-4 500 V，離子源氣體Ⅰ（gasⅠ，GSⅠ）為55 psi，氣體Ⅱ（gasⅡ，GSⅡ）為60 psi，氣簾氣（curtain gas，CUR）為25 psi，碰撞誘導(dǎo)電離（collision-activated dissociation，CAD）參數(shù)設(shè)置為高。檢測結(jié)束后，根據(jù)二級譜信息利用島津株式會社（https://www.shimadzu.com/）提供的化合物比對數(shù)據(jù)庫進(jìn)行花青素定性。以二甲基亞砜（dimethyl sulfoxide，DMSO）為標(biāo)準(zhǔn)品，利用三重四級桿質(zhì)譜的多反應(yīng)監(jiān)測模式（multiple reaction monitoring，MRM）進(jìn)行花青素的定量分析。

1.4 DNA和RNA提取

將保存于-80 ℃的葉片和種皮樣品分別在液氮中充分研磨。每份葉片樣本取0.1 g凍干粉末，用CTAB法提取基因組DNA[38]。每份種皮樣本取0.1 g凍干粉末，用植物總RNA提取試劑盒FastPure?Universal Plant Total RNA Isolation Kit（#RC411, Vazyme）按照產(chǎn)品說明書提取種皮總RNA。用Gene Company Limited（基因有限公司）的NanoDrop-1000分光光度計檢測DNA和RNA的濃度和質(zhì)量，-80 ℃保存?zhèn)溆谩?/p>

1.5 BSA測序分析

1.5.1 文庫構(gòu)建和測序 依據(jù)《大豆種質(zhì)資源描述規(guī)范和數(shù)據(jù)標(biāo)準(zhǔn)》[39]，從RIL群體中選擇30個黃色種皮品系和30個紅色種皮品系，將DNA分別等量混合，構(gòu)成2個子代基因池，命名為Y和R。在親本SN14和TXAJH中各取10個單株，將DNA分別等量混合，構(gòu)成2個親本基因池，命名為SN和TX。

文庫的構(gòu)建和測序由北京百邁克生物科技有限公司完成。試驗流程按照Illumina公司提供的標(biāo)準(zhǔn)protocol執(zhí)行，首先用超聲破碎的方法將DNA隨機打斷成350 bp大小的片段，DNA片段經(jīng)末端修復(fù)、3′端加A、加測序接頭、純化后進(jìn)行PCR擴增，最終通過Illumina HiSeq進(jìn)行測序。親本池測序深度為10×，后代混池測序深度為30×。

1.5.2 變異檢測和關(guān)聯(lián)分析 利用Bcltofastq（v1.8.4）對測序結(jié)果進(jìn)行堿基識別，原始序列過濾得到clean reads，以大豆參考基因組（Wm82.a2.v1）為模板，利用bwa軟件對clean reads進(jìn)行拼接組裝[40]。使用samtools（v1.9）過濾冗余reads。使用GATK的HaplotypeCaller（局部單體型組裝）算法進(jìn)行SNP和InDel的變異檢測[41-42]。使用SnpEff軟件對變異位點集進(jìn)行注釋和預(yù)測[43]。通過歐式距離（euclidean distance，ED）算法和index算法進(jìn)行關(guān)聯(lián)分析[44-45]。ED算法用于計算變異位點與性狀的關(guān)聯(lián)程度，ED值越大表明該位點在兩混池間的差異越大。ED值計算公式如下：

式中，A、C、G、T分別表示堿基A、C、G、T在突變混池中的頻率，A、C、G、T分別表示堿基A、C、G、T在野生型混池中的頻率。

index算法用于尋找混池之間變異位點基因型頻率顯著差異的位點，變異位點與性狀關(guān)聯(lián)度越強，則?index越接近于1。index值計算方法如下：

()=/(+)

()=/(+)

Δ=()-()

式中，（）表示（）池來源于母本的深度，（）表示（）池來源于父本的深度。

對關(guān)聯(lián)分析結(jié)果區(qū)域應(yīng)用BLAST軟件在NR[46]、Swiss-Prot、GO[47]、COG[48]、KEGG[49]等數(shù)據(jù)庫中對候選區(qū)域編碼基因進(jìn)行匹配和注釋。

1.6 基因定位

利用Misa軟件在BSA候選區(qū)域檢索SSR位點，通過SoyBase（https://soybase.org/）網(wǎng)站從大豆參考基因組（Wm82.a2.v1）提取位點上下游200 bp序列，使用NCBI網(wǎng)站（https://www.ncbi.nlm.nih.gov/tools/ primer-blast/）進(jìn)行引物設(shè)計（電子附表1）。將DNA用ddH2O稀釋至50 ng·μL-1后用于PCR反應(yīng)。使用96孔板配置50 μL反應(yīng)體系，包括25 μL 2×Rapid Taq Master Mix（#P222，Vazyme）、1 μL稀釋后的DNA、各2 μL的正向和反向引物（10 μmol·L-1），以及20 μL ddH2O。在BIO-RED T100 PCR儀上進(jìn)行反應(yīng)，反應(yīng)程序為95 ℃ 3 min；95 ℃ 15 s，57 ℃ 15 s，72 ℃ 15 s，35個循環(huán)；72 ℃ 5 min。PCR產(chǎn)物經(jīng)6%聚丙烯酰胺凝膠電泳法分離，銀染后讀取基因型。

1.7 qRT-PCR分析

使用NCBI（https://www.ncbi.nlm.nih.gov/tools/ primer-blast/）網(wǎng)站設(shè)計引物（電子附表2）。每個樣品取1 000 ng總RNA使用逆轉(zhuǎn)錄試劑盒HiScript?III RT SuperMix for qPCR（+gDNA wiper，#R323，Vazyme）進(jìn)行逆轉(zhuǎn)錄，產(chǎn)物經(jīng)無核酸酶水稀釋至1 ng·μL-1后用于qRT-PCR分析。使用Hard-Shell?384孔板（BIO-RAD）配置20 μL反應(yīng)體系10 μL Taq Pro Universal SYBR qPCR Master Mix（#Q712，Vazyme）、5 μL稀釋后的cDNA、各2 μL的正向和反向引物（10 μmol·L-1），以及4.2 μL無核酸酶水。在QuantStudio 7 Flex上進(jìn)行反應(yīng)。反應(yīng)程序為95 ℃ 30 s；95 ℃ 10 s，60 ℃ 30 s，40個循環(huán)；熔解曲線分析為95 ℃ 10 s，60 ℃ 60 s，95 ℃ 15 s。以無模板反應(yīng)作為陰性對照，每個反應(yīng)設(shè)置3個技術(shù)重復(fù)。

1.8 數(shù)據(jù)分析

利用卡方分析檢驗群體的分離比例，計算方法如下：

式中，表示實際頻數(shù)，表示理論頻數(shù)。

使用作為內(nèi)參基因，根據(jù)2-??Ct方法計算qRT-PCR反應(yīng)中目的基因的相對表達(dá)量，采用log2計算以后的FPKM值繪制熱圖。用GraphPad Prism 9（v9.4.1）軟件進(jìn)行統(tǒng)計分析和繪圖。

2 結(jié)果

2.1 種皮花青素含量差異

2.1.1 不同發(fā)育階段的種皮顏色分析 從開花后第30天開始對SN14和TXAJH的種皮顏色進(jìn)行以10 d為間隔的連續(xù)觀察直至種子脫水成熟（第70天）。如圖1所示，在第40天以前，二者種皮均為綠色，從第50天開始，SN14種皮逐漸變黃而TXAJH種皮則出現(xiàn)紅色且逐漸加深，表明導(dǎo)致二者種皮顏色差異的代謝物在第40天以后才開始出現(xiàn)明顯積累。

圖1 SN14和TXAJH不同時期種皮顏色

2.1.2 種皮著色過程花青素含量分析 為了解種皮發(fā)育過程中花青素的含量變化，利用UPLC-ESI- MS/MS方法檢測了SN14和TXAJH在開花后第40、50、60和70天的種皮花青素成分。對SN14和TXAJH各時期種皮樣本的總花青素進(jìn)行層次聚類分析，以了解各組樣本之間花青素的總體差異和組內(nèi)樣本之間的重復(fù)性。如圖2-A所示，花青素總成分的聚類分析將SN14與TXAJH分別聚為一類，TXAJH顯色前后的種皮也分別聚為一類。各組樣本間的重復(fù)性良好，但SN14第60天的1個生物重復(fù)樣本（S_60_2）與第50天的樣本（S_50_1,2,3）聚到一起，表明前者在組間的相似性高于組內(nèi)，因此，將S_60_2樣本剔除再進(jìn)行后續(xù)分析。

在SN14和TXAJH各時期的種皮中共檢測到12種花青素。如圖2-B所示，在種子成熟過程中花青素總含量隨種皮顏色變化而變化。TXAJH種皮內(nèi)花青素含量在著色前較低，第50天出現(xiàn)激增并繼續(xù)升高，在第60天以后含量穩(wěn)定；而SN14種皮內(nèi)花青素含量較低，且隨著種子的成熟持續(xù)下降。這一趨勢與SN14和TXAJH種皮顏色變化情況一致，表明種皮顏色與種皮花青素含量相關(guān)。在被檢測的4個時期中，TXAJH種皮的色素含量始終高于SN14，種皮顯色后這一差異達(dá)到極顯著水平。紅色時期的TXAJH種皮中含量占比最高的3種花青素分別是矢車菊素-3-O-葡萄糖苷（Cyanidin-3-O-glucoside，Cy-3-glu）、芍藥花素-3-O-葡萄糖苷（Peonidin-3-O-glucoside，Pn-3-glu）和牽?；ㄋ?3-O-葡萄糖苷（Petunidin-3-O-glucoside，Pt-3- glu）；而SN14則分別為矢車菊素（Cyanidin，Cy）、天竺葵素-3-O-葡萄糖苷（Pelargonidin-3-O-beta-D- glucoside，Pg-3-glu）和錦葵色素-3-O-葡萄糖苷（Malvidin-3-O-glucoside，Mv-3-glu）。2組花青素的總含量在成熟種子的種皮中相差超過200倍。由此推測Cy-3-glu、Pn-3-glu和Pt-3-glu 3種花青素是導(dǎo)致TXAJH種皮顯示紅色的重要原因。

A：種皮總花青素聚類分析，橫坐標(biāo)為不同樣本之間的距離，縱坐標(biāo)為樣本編號，其中，S：SN14；T：TXAJH，字母后數(shù)字分別表示相應(yīng)時期和生物重復(fù)，不同樣本的組別通過顏色區(qū)分。B：種皮發(fā)育過程花青素含量分析，柱形圖表示花青素含量，餅圖表示70 d時花青素在SN14和TXAJH中的相對占比，同一花青素用相同顏色填充，**表示差異極顯著（P＜0.01）

2.2 BSA-seq初步定位候選區(qū)間

2.2.1 種皮色性狀的遺傳分析 前期以黃色種皮品種SN14和紅色種皮品種TXAJH為親本通過雜交和單粒下傳的方式最終構(gòu)建了F9RIL群體。其中，F(xiàn)1種皮顏色為黃色，表明種皮顏色性狀黃色相對紅色為顯性。在F9RIL群體中，黃色種皮的品系共109個，紅色種皮的品系共79個。如表1所示，卡方值2= 4.556＞2(0.05, 1)=3.841，表明種皮顏色性狀分離比不符合預(yù)期的3﹕1，即紅色種皮基因受2對及以上的等位基因控制。

表1 種皮色性狀遺傳分析

2.2.2 數(shù)據(jù)質(zhì)控和變異位點檢測 對Illumina HiSeq測序得到的原始結(jié)果進(jìn)行數(shù)據(jù)過濾，主要過濾低質(zhì)量和不可信的reads。過濾后共得到110.22 Gb的清潔數(shù)據(jù)（clean read），Q30平均達(dá)到93.26%，樣品與參考基因組平均可比對率為99.22%，基因組平均覆蓋深度為23.5×，基因組堿基被至少覆蓋1次的平均比例為98.44%。以上結(jié)果表明測序隨機性良好，可用于后續(xù)分析。親本及混池樣品測序數(shù)據(jù)統(tǒng)計如表2所示。

根據(jù)Clean Reads的定位結(jié)果，剔除重復(fù)Reads后用GATK進(jìn)行SNP和InDel的變異檢測，每個樣本先各自生成gVCF，再進(jìn)行群體joint-genotype。最后過濾掉低可信度或重復(fù)的位點得到最終的變異位點集，通過對各樣品間的SNP和InDel進(jìn)行Venn統(tǒng)計（圖3），在具有相同種皮顏色的親本和混池中同時存在的位點可能與相應(yīng)種皮顏色的形成有關(guān)。統(tǒng)計結(jié)果顯示，共有2 285 504個SNP變異位點和501 857個InDel變異位點可能與種皮顏色相關(guān)。

SNP_veen：樣品間SNP位點統(tǒng)計結(jié)果；InDel_veen：樣品間InDel位點統(tǒng)計結(jié)果。變異位點的統(tǒng)計只涉及位點的位置而與基因型無關(guān)

2.2.3 關(guān)聯(lián)分析 通過對得到的變異位點集進(jìn)行關(guān)聯(lián)分析以獲得相關(guān)基因的候選區(qū)域（圖4），SNP-ED關(guān)聯(lián)分析得到1個候選區(qū)域，包含1 174個基因；SNP-index關(guān)聯(lián)分析得到3個攜帶基因的候選區(qū)域，共包含1 141個基因；InDel-ED關(guān)聯(lián)分析得到1個攜帶基因的候選區(qū)域，包含1 201個基因；InDel-index關(guān)聯(lián)分析得到1個攜帶基因的候選區(qū)域，包含1 132個基因。所有候選區(qū)域的相關(guān)信息詳見電子附表3。

綜上，對不同關(guān)聯(lián)分析方法得到的SNP和InDel關(guān)聯(lián)區(qū)域結(jié)果取交集（表3），最終得到1個候選區(qū)域，將該區(qū)域作為紅色種皮基因候選區(qū)間。

A：SNP-ED關(guān)聯(lián)結(jié)果；B：SNP-index關(guān)聯(lián)分析的結(jié)果；C：InDel-ED關(guān)聯(lián)結(jié)果；D：InDel-index關(guān)聯(lián)分析的結(jié)果。橫坐標(biāo)為染色體名稱，彩色的點代表每個SNP（或InDel）位點的ED（或Δindex）值，黑色線為擬合后的ED（或Δindex）值，紅色虛線代表顯著性關(guān)聯(lián)閾值

表2 樣品測序數(shù)據(jù)統(tǒng)計

表3 不同關(guān)聯(lián)分析方法得到的共同候選區(qū)域

2.3 標(biāo)記連鎖分析和候選基因預(yù)測

2.3.1 標(biāo)記連鎖分析 為縮小大豆紅色種皮候選基因的范圍，基于BSA候選區(qū)域篩選了27對在親本及子代基因池之間具有多態(tài)性的SSR引物（電子附表1）對RIL群體中188個品系（紅色種皮品系79個，黃色種皮品系109個）進(jìn)行標(biāo)記-表型連鎖分析。共獲得10種單倍型（圖5），其中4種來自于紅色種皮（Red），另6種來自于黃色種皮（Yellow），將候選區(qū)間縮小至標(biāo)記08-223與08-259之間，長度約702 kb。區(qū)間內(nèi)在親本間有37個基因發(fā)生非同義變異。

圖5 種皮顏色候選區(qū)間標(biāo)記連鎖分析

2.3.2 候選基因預(yù)測 在Phytozome（https://phytozome- next.jgi.doe.gov/）網(wǎng)站中對定位區(qū)間內(nèi)的37個基因進(jìn)行功能注釋（表4），、、和在數(shù)據(jù)庫中功能未知，推測可能為假基因或大豆特有新基因。在成功注釋的33個基因中，為MYB轉(zhuǎn)錄因子基因，和為bHLH轉(zhuǎn)錄因子基因，它們可以調(diào)控花青素生物合成途徑相關(guān)基因的表達(dá)；為花青素還原酶1基因，催化花青素向原花青素的轉(zhuǎn)化?；ㄇ嗨厥怯绊懼参锓N皮顏色的重要物質(zhì)，因此，預(yù)測這4個基因可能為大豆紅色種皮候選基因。

表4 候選區(qū)間基因及功能注釋

2.4 基因表達(dá)分析

為探究候選基因在種皮發(fā)育過程中的表達(dá)模式，用SN14和TXAJH在開花后第30、40、50、60和70天的種皮總RNA進(jìn)行qRT-PCR分析（圖6），在種皮發(fā)育的5個時期，候選基因在SN14與TXAJH中的表達(dá)模式存在區(qū)別。在SN14種皮發(fā)育的各時期表達(dá)水平較低，僅在發(fā)育前期出現(xiàn)相對較高的表達(dá)水平，但仍低于同時期的TXAJH；在TXAJH種皮中，各候選基因的表達(dá)水平均較高，且與種皮的花青素含量之間具有相似的變化趨勢。將SN14和TXAJH種皮發(fā)育過程中6種主要花青素的含量與候選基因的表達(dá)水平進(jìn)行關(guān)聯(lián)分析（圖6-B），結(jié)果表明，它們之間存在相關(guān)性。其中，Cy和Mv-3-glu的積累與4個候選基因的表達(dá)水平之間存在較明顯的負(fù)相關(guān)；和的表達(dá)水平與Cy-3-glu、Pt-3-glu、Pg-3-glu和Pn-3-glu的積累之間存在正相關(guān)且具有顯著差異。此外，定位區(qū)間內(nèi)其他33個基因的qRT-PCR結(jié)果顯示，它們在SN14與TXAJH種皮中的表達(dá)模式?jīng)]有明顯區(qū)別，表明相關(guān)基因與種皮花青素積累之間的相關(guān)性較低。

通過查閱文獻(xiàn)[50]選擇了8個花青素生物合成途徑的結(jié)構(gòu)基因進(jìn)行表達(dá)量分析（圖6-C），在絕大多數(shù)時期，TXAJH種皮中這些基因的表達(dá)水平均顯著高于SN14。這表明候選基因、和可能在花青素生物合成途徑中不僅局限于對某一個結(jié)構(gòu)基因產(chǎn)生調(diào)控，而是對通路中多個基因的表達(dá)都具有調(diào)控作用。

3 討論

3.1 紅色種皮的形成與Cy-3-glu、Pn-3-glu和Pg-3- glu的積累有關(guān)

大豆有著悠久的種植歷史，在漫長的馴化與選育過程中形成了豐富多樣的種皮顏色[51]。在多種豆科作物中種皮顏色都被證實與其內(nèi)的多酚類物質(zhì)含量相關(guān)[52]。大豆也不例外，Malen?i?等[53]發(fā)現(xiàn)大豆種皮顏色的深淺與種皮花青素（屬多酚類）的含量有關(guān)。然而，現(xiàn)有大豆種皮花青素的研究多集中于成熟種子，對種子發(fā)育過程中種皮花青素含量的研究較為罕見。本研究通過檢測SN14和TXAJH在4個發(fā)育階段種皮的花青素含量，發(fā)現(xiàn)種皮花青素含量隨著種皮顏色的加深而升高。在種皮發(fā)育過程中，TXAJH的種皮花青素含量始終高于SN14，有趣的是這些差異在綠色種皮階段并不顯著，而從第50天開始變得極顯著。Cy-3-glu、Pn-3-glu和Pg-3-glu是TXAJH種皮花青素中含量占比最高的3種色素，在SN14與TXAJH種皮顏色出現(xiàn)差異后，三者在TXAJH種皮中的含量顯著增加，與此同時它們在SN14種皮中未檢出。結(jié)果表明，TXAJH的紅色種皮可能主要是這三種花青素逐步積累的結(jié)果。

A：候選基因表達(dá)熱圖，采用log2計算以后的FPKM值進(jìn)行作圖，顏色梯度可視化基因相對表達(dá)水平；B：候選基因與種皮主要花青素間的關(guān)聯(lián)熱圖，橫坐標(biāo)為花青素，縱坐標(biāo)為基因，顏色梯度可視化花青素與基因間的相關(guān)性系數(shù)；C：花青素合成途徑相關(guān)基因在種皮各時期的表達(dá)情況。*：＜0.05，**：＜0.01，***：＜0.001

A: Heat map of candidate gene expression, plotted using FPKM values after Log2calculations, with color gradients to visualize relative gene expression levels; B: Heat map of the association between candidate genes and the major anthocyanins of the seed coat, the horizontal coordinate is anthocyanins and the vertical coordinate is genes, the color gradient visualizing the correlation coefficient between the anthocyanins and the genes; C: Expression of genes related to the anthocyanin synthesis pathway at various times in the seed coat. *:＜0.05, **:＜0.01, ***:＜0.001

圖6 候選基因表達(dá)及關(guān)聯(lián)分析

Fig. 6 Candidate gene expression and association analysis

3.2 BSA-seq快速初定位紅色種皮基因候選區(qū)間

在過去較長的一段時間里，基于圖位克隆方法，利用目標(biāo)基因緊密連鎖分子標(biāo)記在染色體上的位置來逐步確定和分離目標(biāo)基因是許多物種確定候選基因的主要方法[54]。但該方法需要在全基因組范圍大規(guī)模篩選標(biāo)記，耗時費力，限制了相關(guān)優(yōu)異基因資源的探索。本研究使用基于全基因組測序的BSA分析方法，通過變異位點與目標(biāo)性狀的關(guān)聯(lián)分析，快速將種皮色基因候選區(qū)域定位于第8染色體（Chr.08：3 530 000—12 280 000 bp），提高了基因定位的效率。此外，測序結(jié)果提供了候選區(qū)域內(nèi)基因序列的變異信息，既可以在進(jìn)一步候選區(qū)間縮小中為標(biāo)記開發(fā)提供依據(jù)，也可以了解變異對相關(guān)基因結(jié)構(gòu)的影響，為候選基因的篩選提供依據(jù)。本研究通過BSA-seq得到的初定位區(qū)間大小為8.6 Mb，包含1 127個基因，區(qū)間長度較大，且在index關(guān)聯(lián)分析結(jié)果中變異位點與性狀的關(guān)聯(lián)強度未超過理論閾值，推測可能是由于構(gòu)建子代極端混池用的樣本量偏小，不能充分消除遺傳背景噪音所致。但這并不影響候選區(qū)域的準(zhǔn)確性，這一現(xiàn)象在其他物種和性狀的研究中曾有先例[55-56]。BSA-seq與標(biāo)記連鎖分析相結(jié)合是解決這一問題行之有效的辦法[57]。Lei等[58]通過BSA-seq和經(jīng)典QTL圖譜將水稻耐鹽QTL從4.17 Mb縮小至222 kb，并鑒定最終候選基因；Zhao等[37]通過BSA-seq將花生紫色種皮基因定位在4.7 Mb區(qū)間，隨后結(jié)合精細(xì)定位將范圍縮小至2.7 kb并最終鑒定到候選基因。因此，BSA-seq結(jié)合標(biāo)記連鎖分析對于定位主要QTL或挖掘目標(biāo)基因是必要的。

3.3 紅色種皮候選基因的預(yù)測與表達(dá)驗證

花青素的生物合成調(diào)控方式在多種植物中得到了廣泛的研究[59-61]。大量研究發(fā)現(xiàn)調(diào)節(jié)植物花青素合成的轉(zhuǎn)錄因子主要為MYB、bHLH和WD40，它們可以單獨發(fā)揮調(diào)控作用，但主要還是以形成MBW復(fù)合物的方式共同進(jìn)行調(diào)控[62]。與用途廣泛的bHLH和WDR轉(zhuǎn)錄因子相比，特異性的MYB轉(zhuǎn)錄因子是決定MBW復(fù)合物靶點的關(guān)鍵因素[63-64]，已在多種植物中分離并鑒定，如草莓[65]、胡蘿卜[66]、海棠[67]和矮牽牛雜交種[68]等。迄今為止，大豆中有4個MYB轉(zhuǎn)錄因子基因（）、（）、（）和（）被鑒定出來與大豆花青素合成相關(guān)[26, 64, 69]。在本研究得到的4個候選基因中，曾被證明為ANR1基因[70]，而、和在此前的種皮顏色研究中未曾被報道。與MYB轉(zhuǎn)錄因子相關(guān)，和與轉(zhuǎn)錄因子bHLH相關(guān)，在種皮著色的5個階段，它們在TXAJH種皮中的表達(dá)水平均高于同時期的SN14，基因表達(dá)與主要花青素含量的關(guān)聯(lián)分析表明二者之間存在極強的相關(guān)性。此外，在種皮著色的各時期，TXAJH種皮中花青素合成途徑相關(guān)結(jié)構(gòu)基因的表達(dá)水平均高于SN14，這一表達(dá)模式與3個轉(zhuǎn)錄因子相關(guān)候選基因相似，暗示著、和可能不止對某一個結(jié)構(gòu)基因產(chǎn)生了調(diào)控作用，而是同時影響了多個基因的表達(dá)。但這一猜測還需更深入的研究證實。

4 結(jié)論

SN14和TXAJH種皮的顏色與種皮中花青素的組成和含量有關(guān)，Cy-3-O-glu、Pn-3-O-glu和Pt-3-O-glu是使TXAJH種皮顯紅色的重要原因。、和在種皮中的表達(dá)水平和花青素的含量之間具有極強的相關(guān)性。、和可能對花青素合成途徑的多個結(jié)構(gòu)基因產(chǎn)生了調(diào)控作用。

[1] Guo B F, Sun L P, Jiang S Q, Ren H L, Sun R J, Wei Z Y, Hong H L, Luan X Y, Wang J, Wang X B, Xu D H, Li W B, Guo C H, Qiu L J. Soybean genetic resources contributing to sustainable protein production. Theoretical and Applied Genetics, 2022, 135(11): 4095-4121.

[2] Fang C, Kong F J. Soybean. Current Biolog, 2022, 32(17): R902-R904.

[3] 宋健. 大豆種皮色相關(guān)基因的圖位克隆及功能解析[D]. 北京: 中國農(nóng)業(yè)科學(xué)院, 2019.

SONG J. Mapping cloning and functional analysis of genes relatead to oybean seed coat color[D]. Beijing: Chinese Academy of Agricultural Sciences, 2019. (in Chinese)

[4] Kovinich N, Saleem A, Rintoul T L, Brown D C W, Arnason J T, Miki B. Coloring genetically modified soybean grains with anthocyanins by suppression of the proanthocyanidin genesand. Transgenic Research, 2012, 21(4): 757-771.

[5] Enaru B, Dre?canu G, Pop T D, Stǎnilǎ A, Diaconeasa Z. Anthocyanins: factors affecting their stability and degradation. Antioxidants, 2021, 10(12): 1967.

[6] Pojer E, Mattivi F, Johnson D, Stockley C S. The case for anthocyanin consumption to promote human health: a review. Comprehensive Reviews in Food Science and Food Safety, 2013, 12(5): 483-508.

[7] Lee Y M, Yoon Y, Yoon H, Park H M, Song S, Yeum K J. Dietary anthocyanins against obesity and inflammation. Nutrients, 2017, 9(10): 1089.

[8] Mattioli R, Francioso A, Mosca L, Silva P. Anthocyanins: a comprehensive review of their chemical properties and health effects on cardiovascular and neurodegenerative diseases. Molecules, 2020, 25(17): 3809.

[9] Liu J Q, Zhou H B, Song L, Yang Z J, Qiu M, Wang J, Shi S L. Anthocyanins: promising natural products with diverse pharmacological activities. Molecules, 2021, 26(13): 3807.

[10] Nomi Y, Iwasaki-Kurashige K, Matsumoto H, TSUDA T, KALT W. Therapeutic effects of anthocyanins for vision and eye health. Molecules, 2019, 24(18): 3311.

[11] Ngamsamer C, Sirivarasai J, Sutjarit N. The benefits of anthocyanins against obesity-induced inflammation. Biomolecules, 2022, 12(6): 852.

[12] de Arruda Nascimento E, de Lima Coutinho L, da Silva C J, de Lima V L A G, Dos Santos Aguiar J.anticancer properties of anthocyanins: a systematic review. Biochimica et Biophysica Acta Reviews on Cancer, 2022, 1877(4): 188748.

[13] Naing A H, Kim C K. Abiotic stress-induced anthocyanins in plants: their role in tolerance to abiotic stresses. Physiologia Plantarum, 2021, 172(3): 1711-1723.

[14] Zhang Q L, Zhai J J, Shao L, Lin W, Peng C L. Accumulation of anthocyanins: an adaptation strategy ofto low temperature in winter. Frontiers in Plant Science, 2019, 10: 1049.

[15] Waseem M, Rong X Y, Li Z G. Dissecting the role of a basic-loop-transcription factor,, under salt and drought stresses in transgenicL.. Frontiers in Plant Science, 2019, 10: 734.

[16] Chunthaburee S, Sakuanrungsirikul S, Wongwarat T, Sanitchon J, Pattanagul W, Theerakulpisut T. Changes in anthocyanin content and expression of anthocyanin synthesis genes in seedlings of black glutinous rice in response to salt stress. Asian Journal of Plant Sciences, 2016, 15(3/4): 56-65.

[17] Ai T N, Naing A H, Yun B W, Lim S H, Kim C K. Overexpression ofenhances anthocyanin accumulation and heavy metal stress tolerance in transgenic. Frontiers in Plant Science, 2018, 9: 1388.

[18] 邱紅梅, 陳亮, 侯云龍, 王新風(fēng), 陳健, 馬曉萍, 崔正果, 張玲, 胡金海, 王躍強, 邱麗娟. 大豆種子顏色遺傳調(diào)控機制研究進(jìn)展. 作物學(xué)報, 2021, 47(12): 2299-2313.

Qiu H M, Chen L, Hou Y L, Wang X F, Chen J, Ma X P, Cui Z G, Zhang L, Hu J H, Wang Y Q, Qiu L J. Research progress on genetic regulatory mechanism of seed color in soybean (). Acta Agronomica Sinica, 2021, 47(12): 2299-2313. (in Chinese)

[19] Wang M, Li W Z, Fang C, Xu F, Liu Y C, Wang Z, Yang R, Zhang M, Liu S L, Lu S J, Lin T, Tang J Y, Wang Y Q, Wang H R, Lin H, Zhu B G, Chen M S, Kong F J, Liu B H, Zeng D L, Jackson S A, Chu C C, Tian Z X. Parallel selection on a dormancy gene during domestication of crops from multiple families. Nature Genetics, 2018, 50(10): 1435-1441.

[20] Cho Y B, Jones S I, Vodkin L O. Nonallelic homologous recombination events responsible for copy number variation within an RNA silencing locus. Plant Direct, 2019, 3(8): e00162.

[21] Todd J J, Vodkin L O. Duplications that suppress and deletions that restore expression from asynthase multigene family. The Plant Cell, 1996, 8(4): 687-699.

[22] Guo Y, Qiu L J. Allele-specific marker development and selection efficiencies for both flavonoid 3’- hydroxylase and flavonoid 3’,5’-hydroxylase genes in soybean subgenus. Theoretical and Applied Genetics, 2013, 126(6): 1445-1455.

[23] Zabala G, Vodkin L. A putative autonomous 20.5 kb-CACTA transposon insertion in anallele identifies a new CACTA transposon subfamily in. BMC Plant Biology, 2008, 8: 124.

[24] Zabala G, Vodkin L O. A rearrangement resulting in small tandem repeats in thegene of white flower genotypes is associated with the soybeanlocus. Crop Science, 2007, 47(S2): s113-s124.

[25] Nitarska D, Stefanini C, Haselmair-Gosch C, Miosic S, Walliser B, Mikulic-Petkovsek M, Regos I, Slatnar A, Debener T, Terefe-Ayana D, Vilperte V, Hadersdorfer J, Stich K, Halbwirth H. The rare orange-red colored Euphorbia pulcherrima cultivar ‘Harvest Orange’ shows a nonsense mutation in a flavonoid 3'-hydroxylase allele expressed in the bracts. BMC Plant Biology, 2018, 18(1): 216.

[26] Vikhorev A V, Strygina K V, Khlestkina E K. Duplicated flavonoid 3’-hydroxylase and flavonoid 3’, 5’- hydroxylase genes in barley genome. PeerJ, 2019, 7: e6266.

[27] Zabala G, Vodkin L O. Methylation affects transposition and splicing of a large CACTA transposon from a MYB transcription factor regulating anthocyanin synthase genes in soybean seed coats. PLoS One, 2014, 9(11): e111959.

[28] Cho Y B, Jones S I, Vodkin L O. Mutations inilluminate epistatic interactions of theandloci leading to saddle seed color patterns in. The Plant Cell, 2017, 29(4): 708-725.

[29] Zou C, Wang P X, Xu Y B. Bulked sample analysis in genetics, genomics and crop improvement. Plant Biotechnology Journal, 2016, 14(10): 1941-1955.

[30] Michelmore R W, Paran I, Kesseli R V. Identification of markers linked to disease-resistance genes by bulked segregant analysis: a rapid method to detect markers in specific genomic regions by using segregating populations. Proceedings of the National Academy of Sciences of the United States of America, 1991, 88(21): 9828-9832.

[31] Zhang G, Angeles E R, Abenes M L P, Khush G S, Huang N. RAPD and RFLP mapping of the bacterial blight resistance gene xa-13 in rice. Theoretical and Applied Genetics, 1996, 93(1): 65-70.

[32] Ma J J, Zhao Y H, Chen H, Fu C, Zhu L, Zhou X M, Xia H, Hou L, Li G H, Zhuang W J, Wang X J, Zhao C Z. Genome-wide development of polymorphic microsatellite markers and their application in peanut breeding program. Electronic Journal of Biotechnology, 2020, 44: 25-32.

[33] Pan J W, Zhou X M, Ahmad N, Zhang K, Tang R H, Zhao H L, Jiang J, Tian M D, Li C S, Li A Q, Zhang X Y, He L Q, Ma J, Li X J, Tian R Z, Ma C L, Pandey M K, Varshney R K, Wang X J, Zhao C Z. BSA?seq and genetic mapping identified candidate genes for branching habit in peanut. Theoretical and Applied Genetics, 2022, 135(12): 4457-4468.

[34] Ramirez-Gonzalez R H, Segovia V, Bird N, Fenwick P, Holdgate S, Berry S, Jack P, Caccamo M, Uauy C. RNA-seq bulked segregant analysis enables the identification of high-resolution genetic markers for breeding in hexaploid wheat. Plant Biotechnology Journal, 2015, 13(5): 613-624.

[35] Zegeye W A, Zhang Y X, Cao L Y, Cheng S H. Whole genome resequencing from bulked populations as a rapid QTL and gene identification method in rice. International journal of molecular sciences, 2018, 19(12): 4000.

[36] Huang P, Jiang H, Zhu C m, Barry K, Jenkins J, Sandor L, Schmutz J, Box M S, Kellogg E A, Brutnell T P. Sparse panicle1 is required for inflorescence development inand maize. Nature Plants, 2017, 3: 17054.

[37] Zhao Y H, Ma J J, Li M, Deng L, Li G H, Xia H, Zhao S Z, Hou L, Li P C, Ma C L, Yuan M, Ren L, Gu J Z, Guo B Z, Zhao C Z, Wang X J. Whole-genome resequencing-based QTL-seq identifiedgene encoding a R2R3-MYB transcription factor controlling peanut purple testa colour. Plant Biotechnology Journal, 2020, 18(1): 96-105.

[38] Kalendar R, Boronnikova S, Sepp?nen M. Isolation and purification of DNA from complicated biological samples. Methods in Molecular Biology, 2021, 2222: 57-67.

[39] 邱麗娟，常汝鎮(zhèn). 大豆種質(zhì)資源描述規(guī)范和數(shù)據(jù)標(biāo)準(zhǔn). 北京: 中國農(nóng)業(yè)出版社, 2006.

Qiu L J, Chang R Z. Descriptorsand data standard for melon (spp.). Beijing: China Agricuture Press, 2006. (in Chinese)

[40] Li H, Durbin R. Fast and accurate short read alignment with burrows-wheeler transform. Bioinformatics, 2009, 25(14): 1754-1760.

[41] McKenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, Garimella K, Altshuler D, Gabriel S, Daly M, DePristo M A. The genome analysis toolkit: a mapreduce framework for analyzing next-generation DNA sequencing data. Genome Research, 2010, 20(9): 1297-1303.

[42] Reumers J, De Rijk P, Zhao H, Liekens A, Smeets D, Cleary J, Van Loo P, Van Den Bossche M, Catthoor K, Sabbe B, Despierre E, Vergote I, Hilbush B, Lambrechts D, Del-Favero J. Optimized filtering reduces the error rate in detecting genomic variants by short-read sequencing. Nature Biotechnology, 2012, 30(1): 61-68.

[43] Cingolani P, Platts A, Wang l L, Coon M, Nguyen T, Wang L, Land S J, Lu X Y, Ruden D M. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome ofstrain w1118; iso-2; iso-3. Fly, 2012, 6(2): 80-92.

[44] Hill J T, Demarest B L, Bisgrove B W, Gorsi B, Su Y C, Yost H J. MMAPPR: mutation mapping analysis pipeline for pooled RNA-seq. Genome Research, 2013, 23(4): 687-697.

[45] Fekih R, Takagi H, Tamiru M, Abe A, Natsume S, Yaegashi H, Sharma S, Sharma S, Kanzaki H, Matsumura H, Saitoh H, Mitsuoka C, Utsushi H, Uemura A, Kanzaki E, Kosugi S, Yoshida K, Cano L, Kamoun S, Terauchi R. MutMap+: genetic mapping and mutant identification without crossing in rice. PLoS One, 2013, 8(7): e68529.

[46] 鄧泱泱, 荔建琦, 吳松鋒, 朱云平, 陳耀文, 賀福初. nr數(shù)據(jù)庫分析及其本地化. 計算機工程, 2006, 32(5): 71-73, 76.

Deng Y Y, Li J Q, Wu S F, Zhu Y P, Chen Y W, He F C. Integrated nr database in protein annotation system and its localization. Computer Engineering, 2006, 32(5): 71-73, 76. (in Chinese)

[47] Ashburner M, Ball C A, Blake J A, Botstein D, Butler H, Cherry J M, Davis A P, Dolinski K, Dwight S S, Eppig J T, Harris M A,Hill D P, Issel-Tarver L, Kasarskis A, Lewis S, Matese J C, Richardson J E, Ringwald M, Rubin G M, Sherlock G. Gene ontology: tool for the unification of biology. The gene ontology consortium. Nature Genetics, 2000, 25(1): 25-29.

[48] Tatusov R L, Galperin M Y, Natale D A, Koonin E V. The COG database: a tool for genome-scale analysis of protein functions and evolution. Nucleic Acids Research, 2000, 28(1): 33-36.

[49] Kanehisa M, Goto S, Kawashima S, Okuno Y, Hattori M. The KEGG resource for deciphering the genome. Nucleic Acids Research, 2004, 32 (Supp_1): D277-D280.

[50] Lu N, Rao X L, Li Y, Jun J H, Dixon R A. Dissecting the transcriptional regulation of proanthocyanidin and anthocyanin biosynthesis in soybean (). Plant Biotechnology Jounal, 2021, 19(7): 1429-1442.

[51] Yuan B Q, Yuan C P, Wang Y M, Liu X D, Qi G X, Wang Y N, Dong L C, Zhao H K, Li Y Q, Dong Y S. Identification of genetic loci conferring seed coat color based on a high-density map in soybean. Frontiers in Plant Science, 2022, 13: 968618.

[52] Amarowicz R, Ronald B P. Legumes as a source of natural antioxidants. European Journal of Lipid Science and Technology, 2008, 110(10): 865-878.

[53] Malen?i? D, Cveji? J, Miladinovi? J. Polyphenol content and antioxidant properties of colored soybean seeds from central Europe. Journal of Medicinal Food, 2012, 15(1): 89-95.

[54] Lindner H, Raissig M T, Sailer C, Shimosato-Asano H, Bruggmann R, Grossniklaus U. SNP-Ratio mapping (SRM): identifying lethal alleles and mutations in complex genetic backgrounds by next-generation sequencing. Genetics, 2012, 191(4): 1381-1386.

[55] Zhao M M, Hu B L, Fan Y W, Ding G M, Yang W L, Chen Y, Chen Y H, Xie J K, Zhang F T. Identification, analysis, and confirmation of seed storability-related loci in dongxiang wild rice (griff.). Genes, 2021, 12(11): 1831.

[56] Liu D, Wei X, Sun D, Yang S, Su H, Wang Z, Zhao Y, Li L, Liang J, Yang L, Zhang X, Yuan Y. An SNP mutation of geneconverts petal color from purple to white in radish (L.). Frontiers in Plant Science, 2021, 12: 643579.

[57] Wambugu P, Ndjiondjop M N, Furtado A, Henry R. Sequencing of bulks of segregants allows dissection of genetic control of amylose content in rice. Plant Biotechnology Journal, 2018, 16(1): 100-110.

[58] Lei L, Zheng H L, Bi Y L, Yang L M, Liu H L, Wang J G, Sun J, Zhao H W, Li X W, Li J M, Lai Y C, Zou D T. Identification of a major QTL and candidate gene analysis of salt tolerance at the bud burst stage in rice (L.) using QTL-seq and RNA-seq. Rice, 2020, 13(1): 55.

[59] Holton T A, Cornish E C. Genetics and biochemistry of anthocyanin biosynthesis. The Plant Cell, 1995, 7(7): 1071-1083.

[60] Winkel-Shirley B. Flavonoid biosynthesis. a colorful model for genetics, biochemistry, cell biology, and biotechnology. Plant Physiology, 2001, 126(2): 485-493.

[61] Jaakola L. New insights into the regulation of anthocyanin biosynthesis in fruits. Trends in Plant Science, 2013, 18(9): 477-483.

[62] Xu W J, Dubos C, Lepiniec L. Transcriptional control of flavonoid biosynthesis by MYB-bHLH-WDR complexes. Trends in Plant Science, 2015, 20(3): 176-185.

[63] Li Y Q, Shan X T, Tong L N, Wei C, Lu K Y, Li S Y, Kimani S, Wang S C, Wang L, Gao X. The conserved and particular roles of theregulatorfromhybrida in flower anthocyanin biosynthesis. Plant and Cell Physiology, 2020, 61(7): 1365-1380.

[64] 高瑞芳.與協(xié)同調(diào)控大豆種皮顏色的機制解析[D]. 長春: 東北師范大學(xué), 2021.

Gao R F. Mechanism analysis ofandsynergistic regulation of soybean seed coat color[D]. Changchun: Northeast Normal University, 2021. (in Chinese)

[65] Zhang Z Y, Shi Y N, Ma Y C, Yang X F, Yin X R, Zhang Y Y, Xiao Y W, Liu W L, Li Y D, Li S J, Liu X F, Grierson D, Allan A C, Jiang G H, Chen K S. The strawberry transcription factor FaRAV1 positively regulates anthocyanin accumulation by activation of FaMYB10 and anthocyanin pathway genes. Plant Biotechnology Journal, 2020, 18(11): 2267-2279.

[66] Xu Z S, Yang Q Q, Feng K, Yu X, Xiong A S. DcMYB113, a root-specific R2R3-MYB, conditions anthocyanin biosynthesis and modification in carrot. Plant Biotechnology Journal, 2020, 18(7): 1585-1597.

[67] Tian J, Peng Z, Zhang J, Song T T, Wan H H, Zhang M L, Yao Y C.regulates coloration via activatingand later structural genes in ever-red leaf crabapple. Plant Biotechnology Journal, 2015, 13(7): 948-961.

[68] Shimada S, Otsuki H, Sakuta M. Transcriptional control of anthocyanin biosynthetic genes in the Caryophyllales. Journal of Experimental Botany, 2007, 58(5): 957-967.

[69] Yan F, Di S, Takahashi R. CACTA-superfamily transposable element is inserted in MYB transcription factor gene of soybean line producing variegated seeds. Genome, 2015, 58(8): 365-374.

[70] Kovinich N, Saleem A, Arnason J T, Miki B. Identification of two anthocyanidin reductase genes and three red-brown soybean accessions with reduced anthocyanidin reductase-1 mRNA, activity, and seed coat proanthocyanidin amounts. Journal of Agricultural and Food Chemistry, 2012, 60(2): 574-584.

Pigment Identification and Gene Mapping in Red Seed Coat of Soybean

CAO Jie1,2, GU YongZhe2, HONG HuiLong2, WU HaiTao2, ZHANG Xia2, SUN JianQiang3, BAO LiGao4, QIU LiJuan1,2

1College of Life Sciences, Jilin Agricultural University, Changchun 130118;2Institute of Crop Sciences, Chinese Academy of Agricultural Sciences, Beijing 100081;3College of Agriculture, Northeast Agricultural University, Harbin 150030;4Agriculture and Animal Husbandry Technology Promotion Center of Inner Mongolia Autonomous Region, Hohhot 010018

【Objective】To identify the key genes controlling anthocyanin synthesis and accumulation, to uncover changes in anthocyanin content of the seed coat during seed development, and the primary anthocyanin components responsible for the red seed coat of Taixingaijiaohong (TXAJH); and to lay the groundwork for a thorough understanding of the regulatory mechanism of red seed coat formation.【Method】Using ultra-high performance liquid chromatography-tandem mass spectrometry (UPLC-ESI-MS/MS), the anthocyanin composition and concentration of the yellow seed coat of soybean Suinong 14 (SN14) and the red seed coat of soybean TXAJH at various developmental stages were identified. The potential areas of red testa-related genes were first identified using bulked segregant analysis (BSA) on the recombinant inbred lines (RILs) made by crossing SN14 and TXAJH. Based on this discovery, we performed marker linkage analysis to restrict the candidate intervals and predict the candidate genes, and qRT-PCR to confirm the expression of the anticipated candidate genes.【Result】When seed coats from the four developmental phases of SN14 and TXAJH were analyzed, a total of 12 anthocyanins were discovered. Cluster analysis of total anthocyanins revealed substantial changes in the seed coat's anthocyanin composition between TXAJH and SN14 as well as between TXAJH before and after color development. The anthocyanin content of the SN14 seed coat gradually decreased as the seed developed, whereas the TXAJH seed coat's content increased quickly and remained stable. After the development of the seed coat's color, the anthocyanin contents of SN14 and TXAJH showed highly significant differences, and at the mature stage, the TXAJH seed coat's anthocyanin content was more than 200 times that of SN14. The crimson coloring of the TXAJH seed coat was largely due to cyanidin-3-O-glucoside (Cy-3-glu), peonidin-3-O-glucoside (Pn-3-glu), and petunidin-3-O-glucoside (Pt-3-glu). The candidate interval for the red seed coat gene on chromosome 8 was discovered at 8.66 Mb by BSA-seq association analysis. 27 polymorphic markers were used in the marker linkage analysis, which produced 10 haplotypes and reduced the candidate interval to 702 kb. Nonsynonymous variations in 37 genes between the parents were found during this interval, these include the genes for encode the anthocyanin reductase 1 (), the bHLH transcription factor (and), and the MYB transcript factor (). These genes may be involved in regulating the biosynthesis of anthocyanins, and anthocyanin reductase 1 can convert anthocyanins to proanthocyanidins (PA). The results of gene expression analysis revealed that candidate genes and genes related to the anthocyanin biosynthesis pathway had comparable expression patterns in SN14 and TXAJH, and both were expressed at lower levels in SN14 and at higher levels in TXAJH. It was discovered that there was a significant link between the principal constituents of seed coat anthocyanins and the level of candidate gene expression.【Conclusion】The anthocyanin makeup of SN14 and TXAJH's seed coats differed, and Cy-3-glu, Pn-3-glu, and Pt-3-glu may be to blame for the TXAJH's seed coat's red hue. According to predictions,,,, andwill likely be a candidate gene for the red seed coat, in which,, andmay control a number of anthocyanin biosynthesis pathway genes.

soybean; seed coat color; anthocyanin; BSA-seq; gene mapping; transcription factors

10.3864/j.issn.0578-1752.2023.14.002

2023-03-03；

2023-04-23

國家重點研發(fā)計劃（2021YFD1201600）、中央級公益性科研院所基本科研業(yè)務(wù)費專項（S2022ZD02）

曹杰，E-mail：cj291@qq.com。通信作者邱麗娟，E-mail：qiulijuan@caas.cn

（責(zé)任編輯 李莉）