王嵐春， 沈方圓， 歐陽丹， 李 校

(四川大學生命科學學院 生物資源與生態(tài)環(huán)境教育部重點實驗室， 成都 610064)

1 引 言

生物體會通過基因的表達和沉默來調控自身的生長發(fā)育，其調控過程是相互獨立而又相互依存的.這種調控可能發(fā)生在基因表達過程中的任意一個階段，如 DNA 的轉錄過程，mRNA的加工過程，以及 mRNA的翻譯過程等，且調控的發(fā)生需要各種酶和調節(jié)蛋白的相互配合[1].其中轉錄因子參與DNA轉錄起始過程來影響基因的表達.轉錄因子通常被定義為具有序列特異性的DNA結合并能夠激活或抑制轉錄的蛋白質[2].基因表達在轉錄水平上的調控影響或控制著細胞或有機體中的許多生物學過程，如細胞周期的進展、代謝和生理平衡以及對環(huán)境的反應.

漫長的自然選擇與進化過程當中，植物形成了其特有的基因表達和調控機制，以適應不斷改變的生存環(huán)境.轉錄因子在這一過程中起到了很大的作用.MYB轉錄因子廣泛分布于所有真核生物中，是植物界最大的轉錄因子家族之一[1].大多數(shù)MYB蛋白起著轉錄因子的作用，具有不同數(shù)量的MYB結構域重復序列，賦予它們結合DNA的能力[1].MYB蛋白是控制發(fā)育、代謝和對生物和非生物脅迫反應的調控網絡中的關鍵因子[3].

MYB轉錄因子(TFs)是一組由1到4個MYB重復序列定義的全真核轉錄因子.50個氨基酸組成三個α螺旋，每個重復序列的第二個和第三個螺旋形成一個螺旋-旋轉-螺旋(HTH)結構，其中有三個間距規(guī)則的色氨酸(或疏水)殘基，形成疏水核心.每個重復序列的第三個螺旋是與DNA直接接觸的DNA識別螺旋[4]. HTH與啟動子中的調節(jié)元件相互作用，而C-末端區(qū)域負責與真核轉錄機制的其他成分建立蛋白質-蛋白質相互作用[5].在DNA接觸過程中，兩個MYB重復序列緊密地堆積在主溝中，使兩個識別螺旋協(xié)同結合到特定的DNA識別序列基序上.MYB-DNA結合域的長度為100～160個殘基，這取決于N-末端區(qū)域不完全重復(稱為R)的數(shù)量.根據(jù)MYB重復數(shù)和MYB重復序列的特性，MYB蛋白質一般分為MYB-related、R2R3-MYB、R1R2R3-MYB和4R-MYB蛋白質[1].在這四種類型中，R2R3MYB是高等植物特有的，在大多數(shù)植物中數(shù)量占優(yōu)勢，其特點是存在保守的MYB結構域和高度可變的C-末端區(qū)域[6, 7].MYB TFs是植物中最大的TF類之一，該類TF的規(guī)模主要歸因于R2R3-MYB TF家族的迅速擴張[8].先前的研究表明， R2R3-MYB TFs的數(shù)量隨著綠色植物進化過程而增加[9].

越來越多的證據(jù)表明，R2R3MYB轉錄因子參與植物許多生理生化過程，如葉毛狀體分化[10]、次生壁形成[11]、花藥和花粉發(fā)育[12, 13]、腋生分生組織形成[14]，二次代謝的調節(jié)包括類黃酮[15]、花青素[16]和木質素[17].除此之外，R2R3MYB家族成員還參與植物對各種非生物和生物脅迫的防御和響應[18-20]，并在調節(jié)植物對包括吲哚乙酸[21]、脫落酸[21, 22]、赤霉素[23, 24]、乙烯、水楊酸和茉莉酸[24]等植物激素線索的反應中發(fā)揮了作用，以及一些環(huán)境信號響應，如水可用性[25]，光[26]和營養(yǎng)元素[27].

隨著多種植物的全基因組測序完成，實驗技術及手段的進步，使得科研人員可以對 R2R3-MYB 基因家族在內的基因資源進行更為廣泛和深入的研究.對于R2R3-MYB轉錄因子的生物功能研究也有了很大的進展，越來越多的R2R3-MYB轉錄因子的功能被揭示.雖然在一些物種中已經對該家族進行了全基因組分析，但對簸箕柳(Salixsuchowensis)的R2R3-MYB基因知之甚少.簸箕柳屬于楊柳科中的柳屬，筐柳組，因其枝條強韌，經常被用來編制筐籃等農具，是一種有發(fā)展前途的經濟植物.除此之外，簸箕柳也可用作固砂樹種、河堤邊的防浪林樹種，具有一定的生態(tài)價值.對簸箕柳進行R2R3-MYB基因的全基因組鑒定，且對其進行該家族的多樣性、進化和其功能研究，有助于了解該蛋白的擴展過程和機制以及隨后該蛋白的亞功能化和新功能化機制.進一步探索R2R3-MYB參與簸箕柳生命生長過程的途徑、方式、調控機制還可以為改良作物遺傳及抗逆提供理論依據(jù)及支持，也可以直接作為基因資源加以利用.更為詳細的關于系統(tǒng)發(fā)育分析、輔助基序和DNA結合特異性的發(fā)現(xiàn)也會為深入了解植物R2R3-MYB轉錄因子的進化史提供線索.除此之外，因簸箕柳枝條強勁、十分耐澇，是河堤邊防浪林樹種的良好選擇，其具有耐性和韌性的枝條在水波中起到了很好的緩沖作用，但洪水或者大風天過后往往造成樹木的倒伏，研究簸箕柳的抗重力刺激響應為后期提升簸箕柳的抗倒伏性能提供了一定的理論基礎.

2 材料與方法

2.1 材 料

從PopGenIE[28]數(shù)據(jù)庫(http://popgenie.org/)下載簸箕柳的基因組、蛋白質序列和注釋文件；從轉錄因子數(shù)據(jù)庫PlantTFDB[29]v5.0(http://planttfdb.gao-lab.org/)下載擬南芥MYB轉錄因子蛋白序列.

2.2 方 法

2.2.1 簸箕柳R2R3-MYB轉錄因子家族的鑒定及其保守結構域分析 根據(jù)Pfam數(shù)據(jù)庫[30](http://pfam.xfam.org/)中MYB轉錄因子家族保守結構域的HMM文件(PF00249)，用HMMER 3.0[31]軟件對簸箕柳蛋白質序列進行本地搜索.將擬南芥R2R3-MYB轉錄因子蛋白序列作為查詢序列，對簸箕柳蛋白質數(shù)據(jù)庫進行BLASTP搜索.將簸箕柳蛋白質序列傳入iTAK數(shù)據(jù)庫[32](http://itak.feilab.net/cgi-bin/itak/index.cgi)進行線上轉錄因子家族鑒定.將上述方法得到的蛋白序列整合、去除冗余后，批量搜索Pfam數(shù)據(jù)庫[30]和NCBI-CDD數(shù)據(jù)庫[33](https://www.ncbi.nlm.nih.gov/Structure/bwrpsb/bwrpsb.cgi)去除結構域不完整的蛋白序列.將鑒定出的蛋白進行重新命名：SsMYB001-SsMYB158.

將上一步鑒定出的簸箕柳R2R3-MYB轉錄因子家族的蛋白序列導入ClustalW[34]進行多序列比對，將得到的結果進行統(tǒng)計.利用網絡服務器MEME[35](http://meme-suite.org/tools/meme)繪制保守結構域sequence logo圖形.

2.2.2 系統(tǒng)發(fā)育、基序組成、基因結構分析 將上一步鑒定出來的簸箕柳R2R3-MYB轉錄因子家族成員的序列導入MAGA-X[36]軟件，采用MAGA-X自帶的ClustalW進行多序列比對，鄰接法(NJ)，檢驗參數(shù)(Bootstrap)設置為1000,進行系統(tǒng)發(fā)育分析,將得到的結果進行亞組分類.使用網絡服務器MEME[35]尋找簸箕柳R2R3-MYB轉錄因子家族成員的保守基序，預期搜索數(shù)量設置為20，允許最小寬度為6，允許最大寬度為50.

使用在線軟件iTOL[37](https://itol.embl.de/)進行系統(tǒng)進化分析和基序組成分析可視化.根據(jù)簸箕柳注釋文件，使用TBtools[38]進行基因結構分析可視化.

2.2.3 亞細胞定位及理化性質分析 使用在線軟件CELLO[39](http://cello.life.nctu.edu.tw/)對簸箕柳R2R3-MYB轉錄因子家族成員進行亞細胞定位預測.采用本地Perl腳本對簸箕柳R2R3-MYB轉錄因子家族成員的基本理化性質進行分析，如理論等電點(pI)、分子量等.

2.2.4 染色體定位和基因復制分析 采用多重共線性工具包MCScanX[40]對簸箕柳R2R3-MYB轉錄因子家族成員進行基因復制事件分析，利用TBtools[38]和Circos[41]工具繪制共線性圖形.采用kaks_Calculator2.0[42]軟件計算基因對的ka/ks值.使用RepeatMasker[43]對簸箕柳基因組進行重復序列鑒定，并統(tǒng)計重復序列和基因密度，基于相對于染色體其他區(qū)域的高重復序列豐度和低基因含量，預測了每條染色體的著絲粒位置.

2.2.5 GO注釋與基因表達分析 將SsMYBs與NCBI的Nr數(shù)據(jù)庫進行blast比對，結果導入Blast2GO[44]，最后使用R語言的GGplot2將注釋結果可視化.

從NCBI_SRA數(shù)據(jù)庫(https://www.ncbi.nlm.nih.gov/sra/)下載簸箕柳重力刺激相關RNA-Seq數(shù)據(jù)(SRA號：SRR9849616、SRR9849619-SRR9849623).運用HISAT2+Stringtie[45, 46]流程計算各基因表達量，再使用DEseq2[47]對兩組數(shù)據(jù)進行差異表達分析，選擇logFC絕對值大于1且FDR小于0.01的基因為差異表達基因.

3 結果分析

3.1 簸箕柳R2R3-MYB轉錄因子家族鑒定及保守結構域分析

在簸箕柳中，經Pfam數(shù)據(jù)庫和NCBI-CDD數(shù)據(jù)庫驗證，去除相關結構域不完整的蛋白序列，一共鑒定出158個簸箕柳R2R3-MYB轉錄因子家族成員.將鑒定出的蛋白進行重新命名：SsMYB001-SsMYB158.

在植物界，R2R3-MYB是MYB蛋白中最大的一個亞群，它含有由兩個相鄰MYB重復序列組成的高度保守的DNA結合域[48].為了探索SsMYBs(R2R3-MYB)成員結構域的特征，提取了SsMYBs成員的MYB結構域，并對其氨基酸序列進行了多重序列比對，結果可以用來檢測SsMYBs蛋白的R2和R3重復序列中每個殘基位置的保守性.通過數(shù)據(jù)統(tǒng)計和MEME可視化得到結構域的一致性序列(圖1)和相關sequence logo(圖2).結果顯示，SsMYBs結構域平均約含有105個氨基酸殘基(圖2)，兩個重復序列中插入或缺失頻率很低.

根據(jù)先前的報道，R2和R3重復序列具有特征性氨基酸，包括一系列分布均勻且高度保守的Trp(W)色氨酸(或疏水氨基酸)殘基[7].在158個

圖1 簸箕柳R2R3-MYB轉錄因子家族的保守結構域橫坐標為重復序列中氨基酸出現(xiàn)的相對位點，縱坐標為該位點占比最大氨基酸的百分比Fig.1 The conserve domain of S.suchowensis R2R3-MYB proteins

圖2 SsMYBs保守結構域sequence logo圖Fig.2 The sequence logo of conserve domain of SsMYBs

SsMYBs蛋白中，R2重復序列含有3個色氨酸殘基，其分別位于第6、26和46位，這些色氨酸殘基形成疏水核心，是植物MYB結構域的典型標志.其中，其中第一個色氨酸殘基有部分被組氨酸(H)絲氨酸等殘基替代，第二個色氨酸殘基有部分被天冬酰胺(N)替代.在R3重復序列中，大多數(shù)成員的第一個色氨酸殘基(位于6)被苯丙氨酸(F)取代，第二個色氨酸殘基(位于25)和第三個色氨酸殘基(位于44)在幾乎所有SsMYBs中都很好地保存，尤其是存在于所有成員中的第二個色氨酸殘基.除了保守的W外，R2重復序列末端還存在十二個高度保守的氨基酸殘基：R-37、G-39、K-40、S-41、C-42等,R3重復序列中的E-10、G-22、A-29等也高度保守(圖2).如圖1所示，SsMYBs結合域中的保守區(qū)域主要位于兩個R重復序列的第二個和第三個色氨酸殘基(每個重復序列中HTH結構域的第三螺旋)之間.SsMYBs結合域中每個R重復序列的第一個和第二個色氨酸殘基(第一螺旋)之間的氨基酸序列相對不保守.除此之外，通過sequence logo(圖2)可以更直觀的看出R3重復片段結構域比R2要更為保守.相比之下，SsMYBs結構域以外的區(qū)域在長度和氨基酸組成方面的保守性較差.

3.2 系統(tǒng)發(fā)育、基序組成、基因結構分析

通過多重序列比對和Bootstrap分析，研究了SsMYBs轉錄因子家族成員的系統(tǒng)發(fā)育分析(圖3).除SsMYB101未能很好地成簇外,其余157個家族成員一共被分為了23個亞組，不同的顏色代表了不同的亞組(圖3).基因結構也可以提供進化信息[49].為了深入了解SsMYBs家族的進化，對158個SsMYBs基因的內含子外顯子的數(shù)目和排列進行了分析.結果表明，同一亞組中的基因通常具有相似的內-外顯子結構(圖3).SsMYBs成員的基因結構存在顯著差異，包括外顯子和內含子的數(shù)目和相對位置.其中，內含子的數(shù)目從0到12不等，大多數(shù)基因有兩個內含子(65%)或一個內含子(23%).在四個SsMYBs基因中發(fā)現(xiàn)了三個內含子：SsMYB069，SsMYB027，SsMYB033，SsMYB053.而在SsMYB059、SsMYB034、SsMYB058和SsMYB046中分別發(fā)現(xiàn)了7、10、11和12個內含子.此外，10個SsMYBs基因被發(fā)現(xiàn)是無內含子的.同一亞群中同源性最高的成員通常具有相同或平行的外顯子/內含子模式，表現(xiàn)出相似的數(shù)量、位置和外顯子長度.例如，C2中的七個SsMYBs擁有兩個外顯子和一個內含子、C4中的兩個SsMYBs沒有內含子、C11中的十個SsMYBs包括三個外顯子和兩個內含子.但是，在C12和C18中各存在一個例外，SsMYB034、SsMYB059的外顯子和內含子的位置以及長度與同組其他成員存在顯著差異，成員間的遺傳相似性較低.值得注意的是，在每個亞組的末端節(jié)點中發(fā)現(xiàn)一對或多對高度同源的SsMYBs，這表明這些蛋白質具有相似的功能.

利用MEME程序預測了SsMYBs蛋白的保守基序，共鑒定出20個不同的基序.這些基序被命名為基序1～20，同一亞組的SsMYB成員具有相似的基序結構.值得注意的是，motif14只在C7亞組中存在，motif8、motif10、motif13只在C2亞組中存在且十分保守，motif12只在C6亞組中存在，motif17只在C10亞組中存在，motif17只在C17亞組中存在，motif15只在C18亞組中存在.這些亞組中特有的motif可能有助于其功能分化.在同一亞組的成員，通常具有相似的基序組成和相似的基因結構(圖3).

3.3 亞細胞定位及理化性質分析

SsMYBs蛋白的長度在164到1041個氨基酸之間，長度普遍位于200～400 aa，平均長度為332 aa.其中，SsMYB099是最長的蛋白質，含有1041個氨基酸；SsMYB152和SsMYB106是最短的蛋白質，各含有164個氨基酸.SsMYBs的分子量在18.3到117.3 kD之間,平均分子量是37.2 kD.SsMYBs的理論等電點在4.61到10.14之間.158個SsMYBs蛋白呈現(xiàn)不規(guī)則的理論等電點特征，75個酸性蛋白(47.5%)和83個堿性蛋白(52.5%).大部分SsMYBs蛋白(87%)的亞細胞定位結果顯示其位于細胞核內的.

圖 3 簸箕柳R2R3-MYB轉錄因子家族的保守結構域的系統(tǒng)發(fā)育樹、motif組成及基因結構

3.4 染色體定位和基因復制分析

片段重復和串聯(lián)復制是導致植物基因家族擴張的兩個主要原因，由同一基因座擴張而來的基因往往具有相似度較高的序列.根據(jù)簸箕柳注釋信息，對SsMYB基因進行染色體定位分析，有7個基因無法定位到簸箕柳19條染色體上，而其余151個基因則不均勻、無規(guī)律地分布于19條染色體上(圖4).其中1號染色體上該家族基因數(shù)目最多，有15個；11號染色體上的最少，只有一條.SsMYB基因在1、2、3號等染色體有較高的密度.相反，一些大的染色體區(qū)域缺乏SsMYB基因，如11號染色體的中央和底部部分、14號染色體的底部.

在SsMYB基因中存在兩組串聯(lián)復制，涉及5個SsMYB基因，分別位于13號染色體和19號染色體.此外，使用MCScanX對SsMYB基因進行片段重復或全基因組重復分析(圖5).共鑒定出83個片段重復對，涉及101個SsMYB基因.在這101個SsMYB基因中，58個基因只出現(xiàn)一次記錄，43個基因存在于不止一次的片段重復事件中.這些結果表明，約(104)65%的SsMYB基因可能是由重復事件產生的，在簸箕柳MYB基因家族的擴展中起主要作用.

Ka/Ks比值用于估計中性突變、純化選擇和有益突變之間的平衡.計算了串聯(lián)和片段復制基因對的Ka/Ks比值.結果表明，所有Ka/Ks比值均小于1，表明簸箕柳MYB復制成員在進化過程中可能經歷了純化選擇壓力.

圖 4 簸箕柳R2R3-MYB轉錄因子家族在染色體上的分布和串聯(lián)復制事件

圖5 簸箕柳R2R3-MYB轉錄因子家族的片段復制事件

3.5 GO注釋

相似的基因在不同物種中，其功能往往保守的.為了探索SsMYBs的功能，對這158個成員進行了GO注釋.GO注釋可以預測一個功能未知基因可能執(zhí)行的分子功能(Molecular Function)、可能處于的細胞組分(Cellular Component)以及可能參與的生物學過程(Biological Process)[44].

簸箕柳MYB成員的GO注釋中分子功能結果顯示(圖6)，大部分成員被預測為起轉錄調控(46)、細胞分化(23)、氣孔運動調節(jié)(7)等功能.對于細胞組分，約96%(151)的成員預測為處于細胞核內.對于生物學過程，大部分成員被預測為參與DNA 結合、轉錄調控區(qū)域特異性DNA結合、DNA結合轉錄因子活性、轉錄協(xié)同調節(jié)活性、序列特異性DNA結合等生物學過程.

圖 6 SsMYB蛋白的GO注釋Fig.6 Gene ontology annotation of SsMYB proteins

3.6 SsMYBs在簸箕柳重力刺激下的差異表達分析

R2R3-MYB家族成員參與植物對各種非生物和生物脅迫的防御和響應[18-20]，并在調節(jié)植物對包括吲哚乙酸[21]、脫落酸[21, 22]、生長素[50]等植物激素線索的反應中發(fā)揮重要了作用.重力是調節(jié)植物生長和發(fā)育的一種普遍輸入，各種植物譜系和器官已經進化不同的機制來調節(jié)相對于重力的生長方向[51].在被子植物中，重力刺激的反應木被稱為張力木，形成于莖的上側，且產生張力，將莖向上拉.張力木是通過維管形成層中細胞分裂速率的增加而產生的，其特征是導水導管元件數(shù)量減少，并且含有凝膠細胞壁層(G層)的特殊張力木纖維被認為是張力產生的核心[52].

圖 7 SsMYB基因在簸箕柳重力刺激下的差異表達分析Fig.7 Expression analysis of SsMYB gene under gravistimulation of S.suchowensis

為了探究SsMYBs在簸箕柳應對重力刺激時可能發(fā)揮的作用，對正常莖和重力刺激下莖的張力木進行差異表達分析，選取logFC>1或<-1和padj<0.01的MYB基因為差異表達基因.共篩選出20個差異表達R2R3-MYB基因(圖7)，上調基因有11個，下調基因9個.其中，SsMYB103為顯著下調基因，其在擬南芥中的同源基因(AtMYB15)為擬南芥激活木質素生物合成基因所必需的調節(jié)因子[53].G層的特征是低木質素、高纖維素[52].推測SsMYB103參與了G層的形成進而推動了簸箕柳莖的背地性響應.

4 討 論

MYB家族是植物最大的轉錄因子家族之一，參與了植物的多種重要生物學過程，如初級和次級代謝、發(fā)育過程、生物和非生物脅迫反應、細胞和器官形態(tài)發(fā)生以及細胞周期控制[3, 54].MYB轉錄因子的進化和功能一直是研究的熱點.MYB轉錄因子已在一些植物物種中被鑒定與分析，如擬南芥、水稻、馬鈴薯和茶樹[55-58],但對于簸箕柳MYB家族的了解卻很少.

本研究在簸箕柳中一共鑒定出158個R2R3-MYB轉錄因子家族成員.對這158個成員先后進行了保守結構域、系統(tǒng)進化、基因結構、基序組成、亞細胞定位、理化性質、染色體定位、基因復制事件以及GO注釋等分析.

多序列比對后，發(fā)現(xiàn)SsMYBs結構域平均約含有105個氨基酸殘基，其中，R2重復序列含有3個色氨酸殘基來形成疏水核心，R3重復序列中，其疏水核心中第一個色氨酸殘基被苯丙氨酸取代，第二個和第三個色氨酸殘基保守性很強.相比之下，SsMYBs結構域以外的區(qū)域在長度和氨基酸組成方面的保守性較差，與之前其他植物該家族結構域的研究一致.隨后進行系統(tǒng)進化分析，除SsMYB101未能很好地成簇外,其余157個家族成員一共被分為了23個亞組.在同一亞組中的基因通常具有相似的內含子-外顯子結構和相似的基序組成.值得注意的是，某些motif只存在于特定的亞組，這些亞組特異性motif可能有助于其功能分化.簸箕柳染色體組裝還不足夠完善，這可能是導致7個SsMYBs無法定位到染色體上的原因，其余151個基因則在19條染色體上不均勻、無規(guī)律地分布.片段重復和串聯(lián)復制是導致植物基因家族擴展的兩個主要原因.在SsMYBs中存在兩組串聯(lián)復制事件,和83個片段重復對，這些結果表明，片段重復和串聯(lián)復制在簸箕柳MYB基因家族的擴展中起主要作用.Ka/Ks比值用于估計中性突變、純化選擇和有益突變之間的平衡.計算了串聯(lián)和片段復制基因對的Ka/Ks比值.結果顯示，所有Ka/Ks比值均小于1，表明簸箕柳MYB復制成員在進化過程中可能經歷了純化選擇壓力.采用Blast2GO[44]對SsMYBs進行GO注釋.大多數(shù)SsMYBs被預測位于細胞核內，約占96%.SsMYBs被預測參與簸箕柳很多的生物學過程，其中46位成員被預測為參與轉錄調控過程、23位成員被預測為參與細胞分化過程等.大多數(shù)SsMYBs被預測起DNA綁定分子功能.除此之外，本研究對簸箕柳莖抗重力刺激應答相關的SsMYB進行了分析.共鑒定出20個差異表達基因.其中，SsMYB103為顯著下調基因，為擬南芥AtMYB15的同源基因，可能通過參與了G層的形成進而推動了簸箕柳莖的背地性響應.

綜上所述，本研究對簸箕柳R2R3-MYB轉錄因子家族進行了全基因組水平的鑒定與分析，為簸箕柳R2R3-MYB基因功能研究和家族特性提供了全面系統(tǒng)的信息，并為改善簸箕柳的生長調節(jié)和抗倒伏性能提供了有價值的信息.