基于RNA-seq的福建柏R2R3-MYB轉錄因子鑒定和分析

周成城, 謝德金, 楊德明, 何天友, 陳凌艷, 鄭郁善,

(1.福建農林大學園林學院；2.福建農林大學林學院, 福建 福州 350002)

轉錄因子(transcription factors, TFs)又稱為反式作用因子，是一類能夠與DNA序列上的啟動子區域直接或間接地發生特異性作用的蛋白，從而激活或抑制目標基因的轉錄[1].轉錄因子家族參與植物生長發育、代謝等各個方面的基因表達調控，具有重要的生物學意義[2].MYB轉錄因子是植物中最早被克隆到的[3]，也是真核生物中功能最多樣的家族[4].根據DNA結合結構域數目的不同，MYB轉錄因子分為1R-MYB、R2R3-MYB、3R-MYB和4R-MYB，其中R2R3-MYB是植物MYB轉錄因子家族中數量最多的一類[5].R2R3-MYB的主要功能是初級和次級代謝物的調節、植物生長發育調節、生物或非生物脅迫響應等[6].如擬南芥AtMYB83和AtMYB46是AtSND1的直接靶點，可誘導觸發AtMYB58、AtMYB63和AtMYB85的表達，這些基因又通過與相關啟動子AC元件相互作用來上調各種木質素合成基因[7-8]；R2R3-MYB的S7、S6和S2亞家族分別對黃酮醇[9]、花青素及原花青素[10]調控；在月季中，RhMYB108轉錄因子作為乙烯和茉莉酸信號接收器，通過激活和增強衰老相關基因表達從而調節花瓣衰老的啟動的模式[11]；擬南芥AtMYB49通過調節角質層形成和抗氧化防御調節擬南芥的耐鹽性[12].R2R3-MYB除了正調節作用，還可以作為負調控因子行使功能.在擬南芥中，AtMYB4能夠抑制C4H的轉錄，從而抑制芥子酸酯的合成[13].以上研究均表明R2R3-MYB轉錄因子對植物的生長發育、物質合成等具有重要的作用.

福建柏[Fokieniahodginsii(Dunn) Henry et Thomas]為柏科福建柏屬，為我國特有種屬，分布在福建、廣東、云南、廣西等省區[14].其樹干通直、樹型優美、四季常青、適應性強、抗性較好、材性穩定、結構細膩，是我國優良的材用-景觀兩用型樹種之一[15].木質素的合成對于木材的質量和樹形有著重要影響，而R2R3-MYB轉錄因子對木質素合成的影響已經有了較為廣泛的研究，如Du et al[16]對157個玉米MYB轉錄因子和125個擬南芥MYB轉錄因子進行分析，結果表明，有4個亞家族的MYB轉錄因子與木質素合成以及次生壁增厚有關，其中過表達AtMYB26能夠使擬南芥木質部增厚；在菊花中，過表達CmMYB8使許多木質素合成組分編碼基因表達下調，植物木質素含量降低[17]，同樣地，銀合歡的LlMYB1在轉基因煙草中過表達也減少了木質素的含量[18].而目前福建柏在轉錄因子調控木質素合成的研究方面處于萌芽階段，本文旨在利用福建柏轉錄組數據，對福建柏R2R3-MYB轉錄因子進行鑒定和生物信息學分析，為后續深入開展R2R3-MYB轉錄因子調控研究提供依據.

1 材料與方法

1.1 試驗材料

試驗材料為35年生生長良好、無病蟲害的福建柏鱗葉、樹皮(包括周皮和韌皮部)和根，于2019年8月采集于福建省福州市永泰縣大湖國有林場的福建柏天然林，使用無菌水把材料清洗干凈并用無菌濾紙擦干后，用錫箔紙包裹立即放進干冰中，帶回實驗室置于-80 ℃冰箱保存備用.3個部位材料在同一區域內選擇年齡相近、長勢相同的3個福建柏，各設置6個生物學重復.

1.2 試驗方法

1.2.1 總RNA的提取和轉錄組測序 福建柏葉片、樹皮和根的總RNA用總RNA試劑盒(天根生物有限公司)進行提取，提取后用瓊脂糖凝膠電泳檢測總RNA的質量和分光光度計(Nanodrop 2000)測定其濃度和純度.并用反轉錄試劑盒(Takara)將總RNA反轉錄合成cDNA第一條鏈，置于-20 ℃冰箱中保存備用.從野外采集回來的樣本立即送往武漢邁特維爾生物有限公司進行轉錄組測序，測序結果供后續分析使用.

1.2.2 福建柏MYB轉錄因子家族的篩選和鑒定 在福建柏轉錄組數據中，提取已有MYB基因家族功能注釋的unigene數據，并提交至CD-HIT Suite在線數據庫 (http://weizhong-lab. ucsd. edu/cdhit_suite/cgi-bin/index.cgi) ，去除冗余的核酸片段, 去冗余后的數據利用在線網站plantTFDB (http://planttfdb.cbi.pku.edu.cn/) 進行MYB轉錄因子篩選與鑒定.

1.2.3 福建柏MYB轉錄因子家族鑒定、理化性質與功能注釋 利用OmicsBox分析工具對所有鑒定出的MYB轉錄因子序列進行氨基酸相似性比對、鑒定和篩選，并對篩選出的MYB轉錄因子進行GO功能分類.篩選出的MYB氨基酸序列N端所含DNA-Binding結構域使用在線預測工具ExPASy-pROSITE(https://prosite.expasy.org/)預測，并根據Dubos et al[4]的方法進行分類，分為1R-MYB(含MYB-related)、R2R3-MYB、3R-MYB、4R-MYB.運用在線處理軟件ExPASy-Protparam (http://web.expasy.org/protparam/)、ExPASy-ProtScale (http://web.expasy.org/protscale/)和NetPho3.1 server (http://www.cbs.dtu.dk/services/NetPhos/)分析預測理化性質和親疏水性.亞細胞定位預測使用在線網站Cell-PLoc(http://www.csbio.sjtu.edu.cn/bioinf/Cell-PLoc-2/)和WoLF PSORT(https://wolfpsort.hgc.jp/)，信號肽分析使用在線網站SignalP(http://www.cbs.dtu.dk/services/SignalP-4.1/)進行分析.

1.2.4 福建柏R2R3-MYB轉錄因子的保守結構域和進化樹分析 利用MEME在線軟件數據庫(http://meme-suite.org/tools/meme)對福建柏R2R3-MYB的R2和R3保守基序進行分析.從在線網站plantTFDB數據庫中下載擬南芥R2R3-MYB轉錄因子氨基酸序列，利用MEGA7.0軟件與福建柏R2R3-MYB轉錄因子共同進行ClustalW序列比對和構建系統進化樹.

1.2.5 福建柏R2R3-MYB轉錄因子的表達模式分析和驗證 根據福建柏轉錄組數據中R2R3-MYB轉錄因子在葉片、樹皮和根3個部位的FPKM(fragments per kilobase of transcript per million fragments mapped)值，利用R語言進行標準化和去中心化并使用heatmap包構建表達熱圖.選取不同表達模式下的R2R3-MYB轉錄因子進行RT-qPCR驗證.利用Primer Premier 5.0軟件進行引物設計(表1)，選擇前期篩選出的Actin7作為內參基因，使用謝德金等[19]的方法進行qPCR擴增，每個樣品做3個平行試驗, 并做3個生物學重復.后用Ct(2-ΔΔCt)法計算福建柏R2R3-MYB轉錄因子在葉片、樹皮、根中的相對表達量，利用Excel 2010和GraphPad 7.0進行統計、繪圖.

表1 引物序列Table 1 Primer sequences

2 結果與分析

2.1 福建柏R2R3-MYB轉錄因子鑒定、理化性質與功能注釋

在線網站plantTFDB初步篩選出145個福建柏MYB轉錄因子，經過氨基酸序列同源比對和保守結構域分析，有122個MYB轉錄因子具有MYB保守結構域且Blast同源性較高(均超過70%).通過分類得到71個R2R3-MYB、43個1R-MYB、6個3R-MYB和2個4RMYB，并依次對其進行編碼FhMYB1～FhMYB122.71個R2R3-MYB的氨基酸數目為100～858個，相對分子質量11 618.52～96 926.56；總平均疏水系數(GRAVY)為-0.1～-1.0，推測都屬于親水蛋白；不穩定系數均高于40，都屬于不穩定蛋白.Cell-PLoc和WoLF PSORT兩個網站的預測結果顯示71個R2R3-MYB轉錄因子均定位于細胞核中.信號肽分析表明全部序列的N端都沒有信號肽剪切位點.

OmicsBox軟件的GO分析結果顯示,71個R2R3-MYB轉錄因子共注釋到364個GO terms, 其中分子功能(molecular function, MF)、生物過程(biological process, BP) 和細胞組成(cellular component, CC) 涉及到的GO terms分別為150、82和132個.其中DNA binding類別中富集個數最多，在所有轉錄因子中均有注釋，說明均含有DNA結合結構域；細胞組分中細胞結構組分(cellular anatomical entity)和胞內組分(intracellular component)均富集到66個R2R3-MYB轉錄因子，生物學過程中細胞學過程(cellular process)也相對較多，富集數目為50個.

2.2 福建柏R2R3-MYB轉錄因子家族保守結構域基序和進化樹分析

R2-MYB和R3-MYB的特征保守基序為：-W-(X19)-W-(X19)-W-和-(F)-(X18)-W-(X18)-W-[4].由圖1可知，在71個福建柏R2R3-MYB轉錄因子的R2-MYB特征中第13位、第33位和第53位氨基酸上出現色氨酸(W, Trp)，都分別間隔19個氨基酸，表明保守基序完全一致；在R3-MYB特征保守基序中，第13位苯丙氨酸 (F, Phe)被異亮氨酸 (I, Ile)取代，在第32位和第51位氨基酸上出現色氨酸，間隔18個氨基酸，符合R3-MYB保守基序特征.進一步分析表明，R2-MYB基序中第3個色氨酸殘基(W, Tyr)前9位氨基酸殘基的比對一致性高，表明該區域的氨基酸序列相對更保守. R3-MYB基序中2個色氨酸殘基中間的18個氨基酸殘基的比對一致性也較高，也表現出相對較高的保守性，如脯氨酸(P, Pro)、甘氨酸(G, Gla)、精氨酸(R, Arg)、蘇氨酸(T, Thr)等.這些保守基序能夠維持福建柏R2R3-MYB轉錄因子的螺旋—轉角—螺旋(H-T-H)結構，進而與DNA更好地結合以行使功能.

R2：R2-MYB；R3：R3-MYB；W：色氨酸；F：苯丙氨酸；X：氨基酸.圖1 福建柏R2R3-MYB轉錄因子的保守基序Fig.1 Conservative motif of the R2R3-MYB of F.hodginsii

基于生物學功能的研究，對擬南芥R2R3-MYB轉錄因子進行了系統的亞分類[4].將福建柏和擬南芥R2R3-MYB轉錄因子進行了多序列比對和構建系統進化樹，并標注上擬南芥R2R3-MYB的亞類注釋(圖2).除亞類S12、S16、S18和S20外，其他亞類均有福建柏R2R3-MYB轉錄因子與擬南芥R2R3-MYB共同聚為一類，其中S5中聚類最多，有9個福建柏R2R3-MYB與AtMYB123聚為一類，可能參與苯丙烷類代謝途徑的調控.有11個福建柏R2R3-MYB轉錄因子和擬南芥AtMYB111聚為一類，AtMYB111為S7分類，但根據擬南芥聚類方式，其未和AtMYB11和AtMYB12聚在一起，因此單獨注釋為W1，S7分類同樣可能參與苯丙烷類代謝途徑的調控.而有14個福建柏R2R3-MYB未和擬南芥R2R3-MYB聚為一類，而單獨聚類為W2亞類.

圖2 福建柏和擬南芥R2R3-MYB進化樹Fig.2 Phylogenetic tree analysis of R2R3-MYB between F.hodginsii and Arabidopsis thaliana

2.3 福建柏R2R3-MYB轉錄因子基因的表達模式分析

對福建柏葉片、樹皮和根的FPKM值進行表達模式熱圖分析(圖3)，結果表明，福建柏R2R3-MYB轉錄因子的表達具有明顯的組織特異性.根據表達模式的差異分為5個類別，其中分組Ⅰ包括FhMYB12、FhMYB23、FhMYB35、FhMYB36等16個轉錄因子，該類別在葉中表達與在樹皮、根中差異顯著，表達趨勢為上調，且在樹皮和葉的表達差異不大；分組Ⅱ包括FhMYB17、FhMYB19和FhMYB54，該類別的表達在根、樹皮、葉中依次呈現階梯式差異，表達趨勢為上調；分組Ⅲ包括FhMYB10、FhMYB22、FhMYB31、FhMYB32等12個轉錄因子，該類別在樹皮中表達與在葉、根中差異顯著，表達趨勢為下調，且在樹皮和葉的表達差異不顯著；分組Ⅳ包括FhMYB7、FhMYB8、FhMYB9、FhMYB11等21個轉錄因子，該類別在根中表達與在葉、樹皮中差異顯著，表達趨勢為下調；分組Ⅴ包括FhMYB1、FhMYB2、FhMYB3、FhMYB4等19個轉錄因子，該類別在葉中表達與在樹皮、根中差異顯著，表達趨勢為下調.

圖3 福建柏R2R3-MYB表達模式Fig.3 Expression pattern of R2R3-MYB in F.hodginsii

為了進一步驗證福建柏R2R3-MYB轉錄因子表達模式的準確性，從5種表達模式中選取具有代表性的轉錄因子(FhMYB8、FhMYB12、FhMYB17、FhMYB20和FhMYB38)，采用RT-qPCR的方法進行驗證(圖4)，結果表明，FhMYB8在根中相對表達量較高，在皮和葉中表達相對較低且兩者差異較小；FhMYB12在根中相對表達量較高，在皮中次之，在葉中表達相對較低；FhMYB17在皮和葉中相對表達量較高，在根中相對表達量較低且與皮、葉中的表達差異較大；FhMYB20在葉中相對表達量較高，在根中接近無表達；FhMYB38在皮中相對表達量較高，3個部位的表達量差異相對較低.5個轉錄因子的RT-qPCR表達結果的總體趨勢與FPKM值的表達模式基本一致，因此驗證了R2R3-MYB轉錄因子的表達模式和組織特異性.

3 討論與結論

R2R3-MYB轉錄因子是MYB轉錄因子家族中的一個大類，對植物的初生、次生代謝物的合成有著重要的作用[20].福建柏不僅材用價值、景觀價值高，其精油的工業和藥用價值也較高，這些特性都離不開福建柏初生、次生代謝物的合成，如木質素、黃酮類、萜類化合物等的合成.為探究福建柏中重要次生代謝物質的合成調控途徑，對福建柏R2R3-MYB轉錄因子進行分析和研究是十分必要的.通過福建柏轉錄組數據庫、plantTFDB數據庫的篩選和比對，鑒定出71個福建柏R2R3-MYB轉錄因子，均含有R2、R3型MYB，且通過結構域分析和序列比對，R2-MYB、R3-MYB中基序中色氨酸殘基前面的氨基酸殘基的比對一致性高，表明該區域的氨基酸序列更加保守和集中.同時也說明 R2R3-MYB轉錄因子N端的保守重復序列能夠更好地形成H(螺旋)-T(轉角)-H(螺旋)結構,使R2R3-MYB轉錄因子蛋白與基因啟動子區結合, 從而控制基因的轉錄水平[21].

在GO功能分類的分析中，預測福建柏R2R3-MYB轉錄因子主要富集于細胞過程和生物過程調控，這符合擬南芥[5]對R2R3-MYB轉錄因子的功能分類.因此進一步和擬南芥R2R3-MYB轉錄因子共同構建系統進化樹，分析可知，除S12、S16、S18和S20外，其他亞類均有福建柏R2R3-MYB轉錄因子與擬南芥R2R3-MYB共同聚為一類，其中S5中聚類最多.在擬南芥中，S5亞家族參與花青素和黃酮類化合物合成途徑[22].其次S13聚類中也較多，擬南芥S13中AtMYB61的研究表明，AtMYB61能夠影響木質素的沉積[23].另外，W1類別中含有較多福建柏R2R3-MYB轉錄因子與AtMYB111聚在一類，AtMYB111原本在擬南芥中和AtMYB11、AtMYB12聚類，但在本研究中沒有聚在S7中，推測可能由于構建進化樹時，未將所有擬南芥R2R3-MYB轉錄因子進行聚類，而是將每個亞家族中的主要轉錄因子進行聚類.S7亞家族在擬南芥中有著非常重要的作用，不僅參與細胞氧化還原過程的調控[24]，也在黃酮類物質的調控中發揮作用[9].W2類別沒有和任何擬南芥R2R3-MYB聚在一類的轉錄因子，該類別與S15亞家族在同一分支中，S15亞家族主要參與植物表皮細胞類型的分化調控[25-26]，W2類別是否符合S15亞家族的功能注釋值得進一步探究.

基因在不同組織中的表達模式和基因的功能是密切相關的[27]，通過研究R2R3-MYB的表達模式，能夠更好地理解R2R3-MYB在調控次生代謝物質合成過程的機理.FPKM表達模式分析表明福建柏R2R3-MYB轉錄因子具有明顯的組織特異性，同時通過RT-qPCR進行了驗證如：FhMYB8和FhMYB20分別在根和葉片中相對表達較高，且進化樹分析都與AtMYB111聚為一類，AtMYB111屬于S7亞家族，推測FhMYB8和FhMYB20可能分別在福建柏根和葉中參與黃酮類物質的合成；FhMYB38在樹皮中表達較高，進化樹分析與S15同一分支，推測FhMYB38可能在福建柏樹皮中調控表皮細胞的分化.

綜上所述，本研究基于福建柏葉、樹皮和根的轉錄組測序數據，篩選和鑒定出71個福建柏R2R3-MYB轉錄因子.對71個R2R3-MYB轉錄因子進行了理化性質、GO富集分類、保守結構域、系統進化樹分析以及表達模式的分析和驗證，能夠初步了解福建柏中R2R3-MYB的數量、功能分類以及表達趨勢，對后續進一步開展福建柏R2R3-MYB轉錄因子對代謝物質的調控具有一定理論價值.