張子璇， 譚明譜， 王 將， 張成才， 羅麗娜， 向增旭，①

(南京農(nóng)業(yè)大學(xué)： a. 園藝學(xué)院， b. 生命科學(xué)學(xué)院， 江蘇 南京 210095)

多倍體植物具有基因劑量效應(yīng)，其體內(nèi)重復(fù)基因的表達(dá)發(fā)生變化，導(dǎo)致植物體在生長發(fā)育、形態(tài)、生理及環(huán)境脅迫耐性等方面表現(xiàn)出優(yōu)勢性狀[1]。多倍體育種是目前藥用植物選育新品種的重要手段，很多藥用植物已成功誘導(dǎo)出多倍體且獲得的多倍體具有活性成分含量高、抗逆性強(qiáng)等特點(diǎn)。例如：穿心蓮〔Andrographispaniculata(Burm. f.) Nees〕同源四倍體葉和莖中穿心蓮內(nèi)酯含量以及莖中脫水穿心蓮內(nèi)脂含量均高于二倍體[2]；甜葉菊(SteviarebaudianaBert.)四倍體的逆境適應(yīng)能力較二倍體更強(qiáng)，能有效利用弱光和低濃度CO2[3]。與二倍體相比，同源四倍體的葉片發(fā)生明顯變化，如PlatycodongrandiflormA. De Candolle同源四倍體的葉片長度和寬度、氣孔長度和寬度、超氧化物歧化酶和過氧化物酶活性、可溶性糖和蛋白質(zhì)含量及總?cè)~綠素含量均升高[4]。

茅蒼術(shù)〔Atractylodeslancea(Thunb.) DC.〕隸屬于菊科(Compositae)蒼術(shù)屬(AtractylodesDC.)，為多年生草本植物，其干燥根莖可入藥，具有燥濕健脾、祛風(fēng)散寒和明目的功效[5]，主要產(chǎn)地為江蘇、湖北及安徽等省份，道地產(chǎn)區(qū)為江蘇鎮(zhèn)江茅山[6]。近年來，由于掠奪性采挖和生境變化等原因，茅蒼術(shù)野生資源瀕臨枯竭，開展茅蒼術(shù)育種和品質(zhì)改良研究迫在眉睫。

轉(zhuǎn)錄組測序技術(shù)能為植物基因表達(dá)研究提供數(shù)據(jù)資料，是研究未完成基因組測序物種的基因及基因功能挖掘的重要手段[7-8]。利用轉(zhuǎn)錄組測序技術(shù)系統(tǒng)分析茅蒼術(shù)二倍體和同源四倍體葉片的轉(zhuǎn)錄組信息，能進(jìn)一步明確茅蒼術(shù)染色體加倍后一系列性狀變異的分子機(jī)制，為挖掘茅蒼術(shù)的功能基因提供基礎(chǔ)資料。本項(xiàng)目組已經(jīng)初步研究了不同倍性茅蒼術(shù)組培苗間的表型差異及生理特征和轉(zhuǎn)錄組差異，結(jié)果表明：與二倍體相比，茅蒼術(shù)同源四倍體葉片發(fā)生明顯的表型變化，包括葉面積指數(shù)、葉長、葉寬、葉厚、葉綠素含量、氣孔的橫徑和縱徑以及保衛(wèi)細(xì)胞中的葉綠體數(shù)顯著增加，葉片的可溶性糖、可溶性蛋白質(zhì)及可溶性淀粉含量升高[9-10]。

鑒于此，以茅蒼術(shù)二倍體和同源四倍體試管苗為研究對象，采用轉(zhuǎn)錄組測序技術(shù)對其葉片進(jìn)行無參轉(zhuǎn)錄組測序分析，并在前期研究[10]基礎(chǔ)上篩選葉片中的高可信度差異表達(dá)unigenes，以揭示茅蒼術(shù)同源四倍體葉片碩大和抗逆性強(qiáng)的分子機(jī)制，并挖掘葉片內(nèi)直接參與同源四倍體葉片性狀變異和抗逆性的相關(guān)候選基因。

1 材料和方法

1.1 材料

供試材料為南京農(nóng)業(yè)大學(xué)園藝學(xué)院中藥生物技術(shù)實(shí)驗(yàn)室前期培養(yǎng)的茅蒼術(shù)二倍體和同源四倍體試管苗。選取實(shí)驗(yàn)室內(nèi)繼代60 d且生長健壯的茅蒼術(shù)二倍體和同源四倍體試管苗，在生根培養(yǎng)基(添加0.5 mg·L-1NAA的1/2MS培養(yǎng)基)中培養(yǎng)20 d后，分別選取長勢基本一致的二倍體和同源四倍體試管苗各5株，每株采集基生葉片2枚，相同倍性的試管苗葉片混勻。按照上述方法重復(fù)取樣3次。采集的葉片經(jīng)液氮速凍后置于-80 ℃冰箱中保存、備用。

1.2 方法

1.2.1 總RNA的提取和純化 采用TRIzol?Reagent(美國Invitrogen公司)分別提取茅蒼術(shù)二倍體和同源四倍體葉片的總RNA；使用NanoDrop 2000超微量分光光度計(jì)(美國Thermo公司)檢測總RNA的濃度和純度；采用DNase Ⅰ〔天根生化科技(北京)有限公司〕純化總RNA，去除總RNA中殘留的基因組DNA；獲得的總RNA經(jīng)Agilent 2100生物分析儀(美國Agilent公司)檢測合格后用于后續(xù)分析。選取質(zhì)量合格的總RNA樣品(質(zhì)量在5 μg及以上，質(zhì)量濃度在200 ng·μL-1及以上，OD260/OD280為1.8～2.2)構(gòu)建cDNA文庫。

1.2.2 cDNA文庫構(gòu)建和轉(zhuǎn)錄組測序 采用磁珠法[11]分離總RNA中的mRNA；采用Fragmentation buffer(美國Invitrogen公司)將獲得的mRNA隨機(jī)斷裂成長度約200 bp的小片段；以mRNA為模板，采用SuperScript Ⅱ Reverse Transcriptase試劑盒(美國Agilent公司)反轉(zhuǎn)錄合成cDNA第1鏈，并使用10×DNA Polymerase buffer和DNA Polymerase Ⅰ試劑(美國Promega公司)合成cDNA第2鏈。使用End Repair Mix(美國Enzymatics公司)將cDNA的粘性末端補(bǔ)成平末端，隨后在3′末端加上1個(gè)A堿基并連接測序接頭；采用Illumina HiSeq 2000高通量測序平臺(tái)(美國Illumina公司)測序，獲得原始測序數(shù)據(jù)。

1.2.3 測序數(shù)據(jù)過濾及組裝 使用Cutadapt v1.16軟件過濾原始測序數(shù)據(jù)，除去污染、含N比高于10%和低質(zhì)量的reads，獲得高質(zhì)量測序數(shù)據(jù)(clean data)。由于茅蒼術(shù)無參考基因組序列，使用Trinity v2.6.6軟件對所有clean data進(jìn)行從頭組裝，去冗余后得到unigenes。

1.2.4 基因功能注釋和差異表達(dá)基因分析 使用BLASTx v2.2.25軟件將得到的unigenes比對到NR、String、SwissProt和KEGG數(shù)據(jù)庫，分別進(jìn)行基因功能注釋；使用Blast2GO軟件(http：∥www.blast2go.com)進(jìn)行GO注釋；使用HMMER v3.2.1軟件進(jìn)行Pfam注釋。

使用edgeR v3.24軟件分析茅蒼術(shù)二倍體和同源四倍體葉片基因的表達(dá)量差異，將log2FC|≥1且FDR≤0.05作為篩選差異表達(dá)基因的標(biāo)準(zhǔn)，其中，F(xiàn)C為差異倍數(shù)，F(xiàn)DR為誤判率。對篩選出的unigenes進(jìn)行GO功能分類統(tǒng)計(jì)和KEGG富集分析；比較項(xiàng)目組前期對茅蒼術(shù)二倍體和同源四倍體試管苗全株[10]及本研究中茅蒼術(shù)二倍體和同源四倍體試管苗葉片的基因差異表達(dá)情況，據(jù)此篩選出葉片中的高可信度差異表達(dá)unigenes。

2 結(jié)果和分析

2.1 轉(zhuǎn)錄組測序數(shù)據(jù)和組裝結(jié)果分析

茅蒼術(shù)二倍體和同源四倍體葉片轉(zhuǎn)錄組的測序分析結(jié)果(表1)表明：原始reads有448 453 428條，clean reads有447 891 390條，所有樣本的Q20均大于98%，Q30均大于94%，GC含量在46%以上。

表1 茅蒼術(shù)二倍體和同源四倍體葉片轉(zhuǎn)錄組測序分析

對獲得的clean reads進(jìn)行組裝和拼接，共獲得750 700個(gè)轉(zhuǎn)錄本數(shù)據(jù)，總堿基數(shù)364 718 682 bp，最長39 024 bp，最短201 bp，平均長度485.8 bp，N50和N90分別為585和239 bp；去冗余后獲得495 797個(gè)unigenes，這些unigenes的總堿基數(shù)為199 717 026 bp，GC含量46.54%，平均長度402.8 bp，N50和N90分別為413和225 bp(表2)。

表2 茅蒼術(shù)葉片轉(zhuǎn)錄本和unigene分析

2.2 基因功能注釋分析

將茅蒼術(shù)轉(zhuǎn)錄組去冗余后獲得的unigenes分別與GO、NR、String、SwissProt、KEGG和Pfam數(shù)據(jù)庫進(jìn)行比對，結(jié)果(表3)顯示：共注釋到259 841個(gè)unigenes，占unigenes總數(shù)的52.41%。其中，注釋到String數(shù)據(jù)庫的unigenes數(shù)最多(180 360)，占unigenes總數(shù)的36.38%；注釋到Pfam數(shù)據(jù)庫的unigenes數(shù)最少(66 311)，占unigenes總數(shù)的13.37%。

表3 茅蒼術(shù)轉(zhuǎn)錄組unigenes與6個(gè)數(shù)據(jù)庫的比對結(jié)果

NR數(shù)據(jù)庫中茅蒼術(shù)與其他植物轉(zhuǎn)錄組unigenes序列的比對結(jié)果(圖1)顯示：茅蒼術(shù)與Cynaracardunculusvar.scolymus(Linn.) Fiori的unigenes序列同源性最高，共有28 614個(gè)unigenes序列相似，占該數(shù)據(jù)庫注釋unigenes數(shù)的20.77%；與向日葵(HelianthusannuusLinn.)的unigenes序列同源性較高，有21 191個(gè)unigenes序列相似，占該數(shù)據(jù)庫注釋unigenes數(shù)的15.38%；與大麥(Hordeumvulgaresubsp.vulgareLinn.)、克里藻〔Klebsormidiumnitens(Kützing) Lokhorst〕、黃石斛(DendrobiumcatenatumLindl.)、旋蒴苣苔〔Boeahygrometrica(Bunge) R. Br.〕、衣藻(Chlamydomonaseustigma)、小立碗蘚〔Physcomitrellapatens(Hedw.) Bruch et Schimp.〕、蓖麻(RicinuscommunisLinn.)和膠球藻C-169(CoccomyxasubellipsoideaC-169)的unigenes序列同源性較低，有2 012～5 475個(gè)unigenes序列相似，均占該數(shù)據(jù)庫注釋unigenes數(shù)的4%以下。另外，茅蒼術(shù)與白梨(PyrusbretschneideriRehd.)、甜菜(Betavulgarissubsp.vulgarisLinn.)和綠藻(Ostreococcustauri)等植物的unigenes序列同源性更低，均占該數(shù)據(jù)庫注釋unigenes數(shù)的1.8%以下。

圖1 NR數(shù)據(jù)庫中茅蒼術(shù)與部分植物轉(zhuǎn)錄組unigenes序列的比對結(jié)果Fig. 1 Comparison results of unigene sequences in transcriptomes of Atractylodes lancea (Thunb.) DC. and some plants in NR database

2.3 差異表達(dá)unigenes分析

比較茅蒼術(shù)二倍體和同源四倍體葉片unigenes的表達(dá)量，篩選出923個(gè)差異表達(dá)unigenes。與二倍體相比，同源四倍體中表達(dá)量上調(diào)的unigenes有317個(gè)，表達(dá)量下調(diào)的unigenes有606個(gè)。

2.3.1 GO功能分類結(jié)果 GO功能分類結(jié)果(表4)表明：茅蒼術(shù)二倍體和同源四倍體葉片差異表達(dá)unigenes的功能被分成生物過程、細(xì)胞組分和分子功能3個(gè)大類61個(gè)小類。在生物過程大類中，注釋為細(xì)胞過程的差異表達(dá)unigenes最多，共有342個(gè)；注釋為代謝過程、應(yīng)激反應(yīng)、生物調(diào)節(jié)和生物過程調(diào)節(jié)的差異表達(dá)unigenes較多，分別有307、249、219和201個(gè)。在細(xì)胞組分大類中，注釋為細(xì)胞、細(xì)胞部分、細(xì)胞器、細(xì)胞器部分和細(xì)胞膜的差異表達(dá)unigenes較多，分別有362、362、337、232和231個(gè)。在分子功能大類中，注釋為結(jié)合和催化活性的差異表達(dá)unigenes較多，分別有291和210個(gè)。

表4 茅蒼術(shù)二倍體和同源四倍體轉(zhuǎn)錄組差異表達(dá)unigenes的GO功能分類

2.3.2 KEGG富集結(jié)果 茅蒼術(shù)二倍體和同源四倍體共有155個(gè)差異表達(dá)unigenes注釋到KEGG通路中，并對其進(jìn)行KEGG分類和通路富集分析，結(jié)果分別見表5和表6。結(jié)果顯示：這些差異表達(dá)unigenes被注釋到代謝、生物系統(tǒng)、遺傳信息處理、環(huán)境信息處理和細(xì)胞過程5個(gè)大類中，以注釋到代謝中的差異表達(dá)unigenes最多(75)，占注釋到KEGG通路的差異表達(dá)unigenes總數(shù)的48.39%。這些差異表達(dá)unigenes富集在8條代謝通路中，包括次生代謝物生物合成(ko01110)、淀粉和蔗糖代謝(ko00500)、脂肪酸延伸(ko00062)、光合作用-天線蛋白(ko00196)、NF-κB信號(hào)通路(ko04064)、植物激素信號(hào)轉(zhuǎn)導(dǎo)(ko04075)、植物MAPK信號(hào)通路(ko04016)和核糖體(ko03010)。

表5 茅蒼術(shù)二倍體和同源四倍體轉(zhuǎn)錄組差異表達(dá)unigenes的KEGG分類

表6 茅蒼術(shù)二倍體和同源四倍體轉(zhuǎn)錄組差異表達(dá)unigenes的KEGG通路富集

2.4 高可信度差異表達(dá)unigenes分析

篩選出的茅蒼術(shù)二倍體和同源四倍體葉片和全株的共有高可信度差異表達(dá)unigenes有298個(gè)，其中，156個(gè)共有高可信度差異表達(dá)unigenes有功能注釋信息，而且其中的58個(gè)在KEGG代謝途徑中有注釋信息，包括熱激蛋白基因(DnaK/Hsp70)、ACC氧化酶基因(ACO)、乙烯受體基因(ETR/ERS)、脂肪酸脫氫酶基因(FAD2)、丙二烯氧化物合酶基因(AOS)、脂氧合酶基因(LOX)、長鏈-3-羥酰基-CoA脫水酶基因(PHS1/PAS2)和14-3-3結(jié)構(gòu)域蛋白基因(YWHAE)。茅蒼術(shù)二倍體和同源四倍體葉片和全株的這些共有高可信度差異表達(dá)unigenes的表達(dá)分析見表7。結(jié)果表明：DnaK/Hsp70、ACO、FAD2和PHS1/PAS2基因的表達(dá)量在同源四倍體中上調(diào)，而ETR/ERS、AOS、LOX和YWHAE基因的表達(dá)量在同源四倍體中下調(diào)。

表7 茅蒼術(shù)二倍體和同源四倍體葉片和全株部分共有高可信度差異表達(dá)unigenes的表達(dá)分析

3 討論和結(jié)論

本研究中，茅蒼術(shù)二倍體和同源四倍體的unigenes在GO、NR、String、SwissProt、KEGG和Pfam數(shù)據(jù)庫的比對結(jié)果表明不是所有unigenes都能在已知數(shù)據(jù)庫中比對出功能，這一結(jié)果與目前測序的許多植物的相關(guān)研究結(jié)果相符[11-13]。這主要是因?yàn)槎鄶?shù)非模式植物尚未完成全基因組測序工作，致使許多植物的基因組信息不全，不可能通過數(shù)據(jù)庫比對得知所有unigenes的功能。值得注意的是，與NR數(shù)據(jù)庫比對結(jié)果顯示，茅蒼術(shù)與同為菊科的C.cardunculusvar.scolymus和向日葵的同源性較高，說明同科植物的unigenes序列同源性更高。

植物多倍化后的性狀差異主要由基因表達(dá)變化引起，且基因表達(dá)只發(fā)生微小改變也可能引起表型變異[14-15]。因此，基因差異表達(dá)是茅蒼術(shù)同源四倍體植株性狀變異的基礎(chǔ)。本研究發(fā)現(xiàn)，注釋為細(xì)胞過程和代謝過程以及細(xì)胞、細(xì)胞部分、細(xì)胞器、細(xì)胞器部分和細(xì)胞膜的茅蒼術(shù)二倍體和同源四倍體葉片的差異表達(dá)unigenes較多，說明茅蒼術(shù)的倍性變化會(huì)影響其細(xì)胞發(fā)育和代謝；還有部分差異表達(dá)unigenes注釋為信號(hào)和信號(hào)傳感器活性，說明茅蒼術(shù)同源四倍體葉片在信號(hào)轉(zhuǎn)導(dǎo)方面與二倍體也存在差異。在其他植物多倍體研究中，酸棗〔Ziziphusjujubavar.spinosa(Bunge) Hu ex H. F. Chow.〕四倍體的差異表達(dá)基因富集在糖類和氨基酸代謝、信號(hào)傳遞以及能量代謝方面，并且這些基因的表達(dá)量在淀粉和蔗糖代謝及植物激素傳導(dǎo)通路中顯著上調(diào)[15]，楊樹(Populusspp.)三倍體中參與糖類代謝和細(xì)胞生長的差異表達(dá)基因的表達(dá)量也上調(diào)[16]。本研究的KEGG富集結(jié)果表明：茅蒼術(shù)二倍體和同源四倍體的差異表達(dá)unigenes富集在脂肪酸延伸、淀粉和蔗糖代謝及次生代謝物生物合成等通路。淀粉、蔗糖和脂肪均為植物細(xì)胞活動(dòng)提供能量，并且，在植物體內(nèi)，淀粉和蔗糖代謝、脂肪酸延伸及次生代謝物生物合成過程均伴隨著能量轉(zhuǎn)化和物質(zhì)積累。分析表明：差異表達(dá)基因可參與植物體的能量轉(zhuǎn)化和物質(zhì)積累過程，為細(xì)胞發(fā)育提供充足的能量和原料、催化相關(guān)酶活性，這也可能是茅蒼術(shù)同源四倍體葉片較二倍體大且厚的一個(gè)重要原因。

植物在適應(yīng)環(huán)境脅迫過程中需要能量、物質(zhì)、激素和酶等共同發(fā)揮作用。在受到環(huán)境脅迫時(shí)，植物體會(huì)產(chǎn)生大量的轉(zhuǎn)導(dǎo)信號(hào)，與淀粉和蔗糖等代謝過程相關(guān)的基因的表達(dá)量明顯上調(diào)，體內(nèi)的淀粉和蔗糖等代謝過程加強(qiáng)[17]，使植物體適應(yīng)逆境條件。在本研究的GO功能分類中，注釋到生物過程類的細(xì)胞過程、代謝過程、應(yīng)激反應(yīng)、生物調(diào)節(jié)和生物過程調(diào)節(jié)等分類的差異表達(dá)unigenes較多，且在KEGG富集結(jié)果中，部分差異表達(dá)unigenes富集于植物激素信號(hào)轉(zhuǎn)導(dǎo)、植物MAPK信號(hào)通路和NF-κB信號(hào)通路共3個(gè)信號(hào)轉(zhuǎn)導(dǎo)相關(guān)通路中。其中，植物MAPK信號(hào)通路是將細(xì)胞膜表面信號(hào)放大并傳遞到細(xì)胞內(nèi)，從而引起細(xì)胞內(nèi)轉(zhuǎn)錄及代謝水平變化的信號(hào)轉(zhuǎn)導(dǎo)途徑之一，與植物生長發(fā)育和抗逆信號(hào)轉(zhuǎn)導(dǎo)具有密切關(guān)系[18]。茅蒼術(shù)同源四倍體葉片在能量轉(zhuǎn)化、物質(zhì)積累和信號(hào)轉(zhuǎn)導(dǎo)方面都與二倍體葉片存在差異，這些差異能夠影響植物體對環(huán)境變化的響應(yīng)，據(jù)此推斷茅蒼術(shù)同源四倍體在應(yīng)對逆境脅迫時(shí)可能存在更強(qiáng)的適應(yīng)能力，以保證整個(gè)植株包括其藥用部位根狀莖的正常生長。后續(xù)研究應(yīng)對不同倍性茅蒼術(shù)藥用部位根莖進(jìn)行轉(zhuǎn)錄組和代謝組聯(lián)合分析，以揭示二者在外觀性狀、生長發(fā)育及抗逆性等方面的分子生物學(xué)差異機(jī)制，為進(jìn)一步定位茅蒼術(shù)的藥用成分合成關(guān)鍵基因及通過基因工程手段進(jìn)行茅蒼術(shù)育種等奠定基礎(chǔ)。

將茅蒼術(shù)二倍體和同源四倍體葉片和全株[10]共有的差異表達(dá)unigenes進(jìn)行比較，共篩選出298個(gè)共有高可信度差異表達(dá)unigenes，其中，葉片高可信度差異表達(dá)unigenes(包括DnaK/Hsp70、ACO、ETR/ERS、FAD2、AOS、LOX、PHS1/PAS2和YWHAE等基因)可能與同源四倍體出現(xiàn)葉變大、變厚，葉色深綠，氣孔增大等表型變化和植株抗逆性相關(guān)[19-34]。本研究中，DnaK/Hsp70、ACO、FAD2和PHS1/PAS2基因的表達(dá)量在茅蒼術(shù)同源四倍體葉片和全株中均上調(diào)，而ETR/ERS、AOS、LOX和YWHAE基因的表達(dá)量在茅蒼術(shù)同源四倍體葉片和全株中均下調(diào)，這可能是茅蒼術(shù)同源四倍體植株抗逆性強(qiáng)于二倍體的原因。另外，茅蒼術(shù)二倍體和同源四倍體葉片和全株[10]非共有的差異表達(dá)unigenes也值得研究，這些基因可能調(diào)控其植株葉片和地下部(如根狀莖)的生長發(fā)育。

本研究結(jié)果顯示：茅蒼術(shù)同源四倍體葉片在細(xì)胞發(fā)育、代謝、能量轉(zhuǎn)化、物質(zhì)積累及信號(hào)轉(zhuǎn)導(dǎo)方面與二倍體間存在差異，差異表達(dá)unigenes可能參與葉片發(fā)育調(diào)控。高可信度差異表達(dá)基因DnaK/Hsp70、ACO、FAD2和PHS1/PAS2的表達(dá)量在茅蒼術(shù)同源四倍體葉片和全株中均上調(diào)，而ETR/ERS、AOS、LOX和YWHAE的表達(dá)量均下調(diào)，這些基因直接參與茅蒼術(shù)同源四倍體葉片發(fā)育和抗逆性，可作為茅蒼術(shù)倍性育種的關(guān)鍵候選基因。