鎘脅迫下菜心的轉(zhuǎn)錄組分析

張尉欣，劉 陽，劉 帥，唐玉林，2

1)深圳大學(xué)生命與海洋科學(xué)學(xué)院，廣東省植物表觀遺傳學(xué)重點(diǎn)實(shí)驗(yàn)室，廣東深圳 518060；2)深圳市海洋生物資源與生態(tài)環(huán)境重點(diǎn)實(shí)驗(yàn)室，廣東深圳 518060

鎘(cadmium，Cd)因其生物有效性高，與其他重金屬(Cu、Pb和As)相比更容易被植物吸收，2001年被聯(lián)合國環(huán)境規(guī)劃署列為全球性危害化學(xué)物質(zhì)之首．Cd脅迫下植株會(huì)出現(xiàn)矮化、發(fā)育推遲和產(chǎn)量減少等毒害效應(yīng)[1]．植物通過調(diào)控許多生物學(xué)過程而產(chǎn)生脅迫響應(yīng)，如金屬硫蛋白和植物螯合肽的螯合作用、Cd2+的區(qū)室化、硫的同化作用、抗氧化機(jī)制的激活、細(xì)胞壁組分對Cd的固定或胞外分泌物與Cd2+的結(jié)合等[2-3]．這些過程的發(fā)生離不開相關(guān)基因的表達(dá)調(diào)控，植物響應(yīng)Cd脅迫的機(jī)理研究已在擬南芥、水稻等植物中有報(bào)道[4-5]，但調(diào)控網(wǎng)絡(luò)有待深入認(rèn)識．

菜心，又名菜薹(Brassicarapassp.chinensisvar.parachinensis)，為十字花科蕓苔屬植物，是一種對鎘不敏感且高積累的植物[6]，也是中國廣東地區(qū)主要種植的蔬菜，而南方地區(qū)土壤的鎘污染較為嚴(yán)重，這勢必造成菜心的鎘污染風(fēng)險(xiǎn)．然而，關(guān)于菜心應(yīng)對Cd脅迫及Cd的吸收、轉(zhuǎn)運(yùn)及積累機(jī)制目前尚不明確．本研究以菜心為研究對象，通過對菜心Cd2+處理前后的葉和根部進(jìn)行轉(zhuǎn)錄組測序和差異表達(dá)基因分析，獲得Cd脅迫后的差異表達(dá)基因，為研究菜心的鎘脅迫響應(yīng)機(jī)制，篩選、培育低Cd積累的菜心品種提供理論依據(jù)．

1 材料與方法

1.1 試供材料

以萌發(fā)后水培15 d的菜心品種“油青49”幼苗為材料，以50 μmol/L Cd(NO3)2的霍格蘭營養(yǎng)液處理幼苗．以處理前(CK)的材料為對照組，處理48 h的材料為處理組，每3株菜心為1個(gè)生物學(xué)重復(fù)，設(shè)3個(gè)生物學(xué)重復(fù)，分別取葉部(L1、L2、L3)和根部(R1、R2、R3)，液氮速凍后，置于-80 ℃冰箱保存．

1.2 cDNA文庫構(gòu)建

用Trizol提取樣品總核糖核苷酸(ribonucleic acid，RNA)，構(gòu)建信使RNA(messenger RNA，mRNA)測序文庫并送北京諾禾致源生物信息科技有限公司進(jìn)行測序，獲得轉(zhuǎn)錄組數(shù)據(jù)．

1.3 轉(zhuǎn)錄組生物信息學(xué)分析流程

對轉(zhuǎn)錄組原始數(shù)據(jù)進(jìn)行生物信息學(xué)分析的流程如下：① 數(shù)據(jù)質(zhì)量控制；② 與參考基因組比對映射；③ 計(jì)算基因表達(dá)量；④ 差異表達(dá)分析；⑤ 差異表達(dá)基因的富集分析．具體方法如下．

1.3.1 轉(zhuǎn)錄組數(shù)據(jù)與參考基因組的比對

原始數(shù)據(jù)利用Trim Galore(v 0.4.1)進(jìn)行質(zhì)控分析，清除接頭和測序質(zhì)量較低的reads等，得到clean reads后使用HISAT2[7](v 2.0.5)與參考基因組大白菜Brassicarapa(IVFCAASv1)[8]進(jìn)行比對映射，比對結(jié)果使用gnm_region_dist(v2.8)計(jì)算讀段(reads)并使用TBtools[9]量化基因表達(dá)量，以FPKM (fragments per kilobase of exon per million fragments mapped)表示．

1.3.2 差異表達(dá)基因的篩選

使用軟件edgeR[10]進(jìn)行基因的差異表達(dá)分析，分別將葉部和根部的處理組和對照組相比較，獲得兩個(gè)差異表達(dá)基因(differential expressed gene，DEG)的集合．采用Benjamini Hochberg校正方法對P值進(jìn)行校正得到校正后的P值，以差異表達(dá)倍數(shù)(fold change， FC)和錯(cuò)誤發(fā)現(xiàn)率(false discovery rate， FDR)作為差異表達(dá)基因的篩選條件，篩選閾值為|lb FC|≥1，F(xiàn)DR<0.05．

1.3.3 差異表達(dá)基因富集分析

將全部基因的蛋白序列輸入到AgBase v 2.0(http://agbase.arizona.edu/index.html)進(jìn)行GO注釋，同時(shí)利用KOBAS 3.0(http://kobas.cbi.pku.edu.cn/index.php)進(jìn)行KEGG注釋，注釋結(jié)果整理后用于后續(xù)富集分析．使用WEGO(http://wego.genomics.org.cn/)和TBtools分別進(jìn)行差異表達(dá)基因的KEGG和GO富集分析．

1.3.4 實(shí)時(shí)熒光定量PCR

從轉(zhuǎn)錄組結(jié)果中選擇9個(gè)差異表達(dá)基因進(jìn)行實(shí)時(shí)熒光定量PCR(real-time quantitative PCR，qRT-PCR)，使用的儀器以及試劑分別為：StepOnePlusTM(Thermon Fisher)、SYBR premix ex taq kit(Takara)．qRT-PCR驗(yàn)證引物序列如表1所示(內(nèi)參基因actin引物序列參考文獻(xiàn)[11])，將轉(zhuǎn)錄組數(shù)據(jù)與qRT-PCR結(jié)果使用軟件Excel進(jìn)行相關(guān)性分析，計(jì)算相關(guān)性系數(shù)R，驗(yàn)證轉(zhuǎn)錄組數(shù)據(jù)的可靠性．

表1 qRT-PCR驗(yàn)證所用引物序列

2 結(jié)果與分析

2.1 轉(zhuǎn)錄組數(shù)據(jù)與參考基因組的序列比對

將對照組葉部(CK_L)和根部(CK_R)， Cd2+處理48 h的葉部(48H_L)和根部(48H_R)進(jìn)行轉(zhuǎn)錄組高通量測序．共獲得原始數(shù)據(jù)63.12 Gbit，去除質(zhì)量低的reads后，各個(gè)樣品過濾后得到的數(shù)據(jù)量達(dá)到3.7 Gbit以上，Q20值(堿基質(zhì)量值高于20的百分比)均大于96%(表2)，說明測序數(shù)據(jù)質(zhì)量較好，能進(jìn)行后續(xù)分析．

利用與菜心同種但不同亞種的大白菜(Brassicarapasubsp.pekinensis)基因組(IVFCAASv1)[8]作為參考基因組進(jìn)行測序片段的比對．轉(zhuǎn)錄組數(shù)據(jù)與該參考基因組的比對率為88.5%～93.2%，唯一比對率為80.53%～89.87%．測序片段能夠較好地比對到參考基因組上，表明測序結(jié)果可進(jìn)行后續(xù)的差異表達(dá)分析．

表2 轉(zhuǎn)錄組數(shù)據(jù)與大白菜參考基因組比對分析

2.2 差異表達(dá)基因的數(shù)量及分布

利用edgeR[10]進(jìn)行對照組和Cd處理組的差異表達(dá)分析，分別獲得葉部和根部的差異表達(dá)基因集合．共3 854個(gè)差異表達(dá)基因，其中，Cd處理后葉部表達(dá)上調(diào)的基因有951個(gè)，下調(diào)的基因有1 269個(gè)；Cd處理后根部表達(dá)上調(diào)的基因有1 086個(gè)，下調(diào)的基因有754個(gè)．

2.3 差異表達(dá)基因的GO分類和富集分析

所有差異表達(dá)基因(3 854個(gè))中，共有2 925個(gè)基因有GO注釋信息．圖1使用WEGO對2 925個(gè)差異表達(dá)基因進(jìn)行GO分類，可見在細(xì)胞組分、細(xì)胞和細(xì)胞部分所占比例較高，其次是膜和細(xì)胞器．在分子功能中，結(jié)合和催化活性所占比例較高，其次是翻譯調(diào)節(jié)活性和轉(zhuǎn)運(yùn)活性．在生物學(xué)過程中，細(xì)胞過程和代謝過程所占比例較高，其次是生物調(diào)節(jié)．

圖1 差異表達(dá)基因GO注釋以及分類分析Fig.1 GO annotation and classification of DEGs

進(jìn)一步對所有差異表達(dá)基因進(jìn)行GO富集分析.選取校正后的P<0.05的富集結(jié)果并按富集系數(shù)進(jìn)行排序，并對排序前20的GO條目進(jìn)行分析，富集選中的基因數(shù)為M， 背景基因數(shù)為N， 結(jié)果見表3．差異表達(dá)基因主要富集在與次生代謝物合成(GO號為0018160、0019419和0006553)、氧化還原反應(yīng)(GO號為0004324、0016661和0004506)、光合作用(GO號為0009538、0009654和0009773)、鐵離子穩(wěn)態(tài)(GO號為0006879和0055072)、木質(zhì)素或木栓質(zhì)的合成及代謝(GO號為0010345和0046274)等有關(guān)的GO條目中．

木質(zhì)素和木栓質(zhì)是細(xì)胞壁中的重要組分，會(huì)影響Cd2+的吸收和運(yùn)輸[12-13]，而Cd2+與Fe2+則可能通過相同的轉(zhuǎn)運(yùn)蛋白進(jìn)行吸收和運(yùn)輸，并與細(xì)胞內(nèi)鐵離子的穩(wěn)態(tài)有關(guān)[14]．進(jìn)一步對富集在木質(zhì)素水解、木栓質(zhì)的生物合成、鐵離子穩(wěn)態(tài)中的相關(guān)基因進(jìn)行分析(圖2)，在木質(zhì)素水解過程中富集到的差異表達(dá)基因主要是漆酶(laccase，LAC)，除了LAC5、LAC7以及LAC12基因以外，其他大部分富集到的LAC基因在Cd脅迫后表達(dá)下調(diào)；而與木栓質(zhì)的生物合成過程相關(guān)的CYP86、KCS2、FAR4、ABC基因在根部表達(dá)上調(diào)，在葉部表達(dá)下調(diào)或無明顯變化；與維持細(xì)胞中鐵離子穩(wěn)態(tài)相關(guān)的基因如FER3和FER1基因在葉部和根部均表達(dá)下調(diào)，而ABC1、BHLH038和Nramp4基因則表達(dá)上調(diào)．

2.4 差異表達(dá)基因KEGG富集分析

對根和葉部中差異表達(dá)均上調(diào)和均下調(diào)基因分別進(jìn)行KEGG通路的富集分析(P<0.05)， 見表4和表5．差異表達(dá)上調(diào)的基因主要富集在植物病原菌互作、MAPK信號傳導(dǎo)通路和植物激素信號傳導(dǎo)、木栓質(zhì)、黃酮以及萜烷等物質(zhì)的生物合成、谷胱甘肽以及硫等物質(zhì)代謝等通路中．差異表達(dá)下調(diào)的基因主要富集在光合作用、核糖體、生物素代謝、脂肪酸的生物合成、芥子油苷的生物合成等通路中．這些代謝通路的變化反映了受Cd脅迫的影響，菜心中與基礎(chǔ)代謝相關(guān)的通路普遍下調(diào)．另一方面，與抵抗脅迫相關(guān)的信號轉(zhuǎn)導(dǎo)通路的激活、次生代謝物的合成與代謝改變，亦是菜心應(yīng)對Cd脅迫的一種主動(dòng)防御方式．

表3 差異表達(dá)基因GO富集分析

1)為富集選中的基因數(shù)/背景基因數(shù)

圖2 GO富集分析中部分GO亞分類的相關(guān)基因的表達(dá)趨勢Fig.2 Gene expression trends of partial GO terms in GO significant enrichment analysis

通 路富集程度校正后P值植物病原菌互作/ plant-pathogen interaction2.101 1.22×10-7MAPK信號通路-植物/ MAPK signaling pathway-plant2.749 8.60×10-7角質(zhì)、木栓質(zhì)和蠟質(zhì)生物合成/ cutin, suberine and wax biosynthesis4.265 5.45×10-5谷胱甘肽代謝/ glutathione metabolism2.969 5.21×10-5苯丙烷類生物合成/ phenylpropanoid biosynthesis2.138 4.34×10-4硫代謝/sulfur metabolism3.681 4.34×10-4甘油脂代謝/ glycerolipid metabolism2.629 0.002 甘油磷脂代謝/ glycerophospholipid metabo-lism2.329 0.004 黃酮和黃酮醇的生物合成/ flavone and flavonol biosyn-thesis6.702 0.005 倍半萜類化合物和三萜類化合 物的生物合成/sesquiterpe-noid and triterpenoid biosynthesis4.160 0.018 玉米素生物合成/ zeatin biosynthesis3.351 0.027 植物激素信號轉(zhuǎn)導(dǎo)/plant hormone signal transduction1.484 0.037 半胱氨酸和蛋氨酸代謝/cysteine and methionine metabolism1.942 0.038

表5 差異表達(dá)下調(diào)基因KEGG富集分析

2.5 差異表達(dá)基因qRT-PCR結(jié)果分析的相關(guān)性

為驗(yàn)證轉(zhuǎn)錄組分析結(jié)果的可靠性，從差異表達(dá)基因中選取9個(gè)與離子轉(zhuǎn)運(yùn)相關(guān)的基因(YSL2、HMA2、ABCG40、CAX1、IRT1、Nramp1、MRP7、ZIP2和NRAMP4)進(jìn)行qRT-PCR分析．結(jié)果顯示，轉(zhuǎn)錄組數(shù)據(jù)與qRT-PCR結(jié)果基本一致，R=0.802 2(圖3)，說明轉(zhuǎn)錄組分析結(jié)果可靠．

圖3 轉(zhuǎn)錄組數(shù)據(jù)與qPCR驗(yàn)證的相關(guān)性分析Fig.3 Correlation analysis of qPCR and transcriptome (RNA-seq) results

3 討 論

經(jīng)檢索，至今尚缺菜心全基因組測序的研究報(bào)道，但進(jìn)化研究均表明其屬于蕓薹(Brassicarapa)的一個(gè)亞種[12]．目前，蕓薹已經(jīng)以大白菜Chiifu-401-42(B.rapassp.pekinensisline Chiifu-401-42)為材料完成了全基因組測序，并繪制了基因組草圖[8]．本研究以該基因組為參考基因組，將菜心的測序數(shù)據(jù)比對到該基因組上，比對效率較高，利用比對結(jié)果進(jìn)行差異表達(dá)分析，對挖掘到的差異表達(dá)基因進(jìn)行qRT-PCR驗(yàn)證，說明利用大白菜的基因組作為參考基因組進(jìn)行菜心測序數(shù)據(jù)分析可行．

Cd脅迫下菜心的轉(zhuǎn)錄組分析結(jié)果表明，在Cd脅迫下差異表達(dá)基因主要富集到的GO條目包括：代謝過程、氧化還原反應(yīng)、光合作用和鐵離子穩(wěn)態(tài)等，而KEGG主要富集的通路為光合作用、谷胱甘肽的代謝、硫代謝、核糖體和脂肪酸的生物合成等．這與水稻[16]、玉米[17]、杞柳[18]和白菜[19]等植物的Cd脅迫下轉(zhuǎn)錄組研究結(jié)果基本一致．

菜心中木栓質(zhì)生物合成、木質(zhì)素代謝和鐵離子穩(wěn)態(tài)的相關(guān)基因的表達(dá)也明顯受到Cd脅迫影響．木栓質(zhì)生物合成途徑中的差異表達(dá)基因有KCS2、FAR4， 及兩個(gè)ABC轉(zhuǎn)運(yùn)家族蛋白基因，它們分別與擬南芥中的AtKCS2、AtFAR4、AtABCG6和AtABCG2基因同源．?dāng)M南芥的AtKCS2[20]、AtFAR4[12]、AtABCG6以及AtABCG2[21]基因均與木栓質(zhì)合成有關(guān)．推測Cd脅迫后，這些基因在根部的上調(diào)表達(dá)會(huì)使木栓質(zhì)合成和含量增加，而木栓質(zhì)是細(xì)胞壁中的具保護(hù)功能的疏水屏障[12]，這一屏障的增厚可能會(huì)減少Cd2+進(jìn)入菜心根細(xì)胞的量，減少對細(xì)胞的毒害．這些基因的差異表達(dá)意味著菜心根部在受到Cd脅迫后啟動(dòng)了這種保護(hù)機(jī)制．

在差異表達(dá)富集的木質(zhì)素代謝途徑的條目中，LAC11、LAC17和LAC4(IRX12)等基因已被證實(shí)參與細(xì)胞壁木質(zhì)素多聚物合成[22-24]．推測Cd脅迫后這些基因的下調(diào)可能會(huì)影響木質(zhì)素的沉積，木質(zhì)素的沉積減弱又會(huì)影響根部次生內(nèi)皮層的發(fā)育，進(jìn)而使根吸收的Cd2+向地上部分的轉(zhuǎn)運(yùn)更順暢[13]，這些基因表達(dá)變化是導(dǎo)致菜心植株在Cd脅迫后Cd含量普遍高于其他十字花科植物的可能原因之一．由此推測，Cd脅迫下菜心可能通過對細(xì)胞壁物質(zhì)(木栓質(zhì)和木質(zhì)素等)的合成與代謝的協(xié)同調(diào)控，提高了其Cd耐受和積累能力，需進(jìn)一步對基因功能和機(jī)理進(jìn)行深入探討．

已有研究表明，Cd會(huì)破壞植物的必需微量元素如鐵、鋅和銅等元素的動(dòng)態(tài)平衡[25]．圖2表明，在Cd脅迫下，菜心中一些與鐵離子穩(wěn)態(tài)相關(guān)的基因的表達(dá)與擬南芥缺鐵時(shí)的表達(dá)趨勢一致：在缺鐵時(shí)，AtNramp4基因表達(dá)上調(diào)[14]，F(xiàn)ER基因表達(dá)下調(diào)[26]，bHLH038基因表達(dá)上調(diào)[27]．也有研究表明，AtNramp3和AtNramp4蛋白除轉(zhuǎn)運(yùn)Fe2+外，還對Cd2+具有極高的轉(zhuǎn)運(yùn)能力[14]．由此推測，Cd2+脅迫時(shí)菜心對Cd2+的吸收影響了細(xì)胞對鐵離子的吸收及細(xì)胞內(nèi)鐵離子穩(wěn)態(tài)，造成植物出現(xiàn)了類似缺鐵的響應(yīng)，而轉(zhuǎn)運(yùn)蛋白Nramp4基因等在葉部表達(dá)上調(diào)則有可能促進(jìn)葉部Cd的轉(zhuǎn)運(yùn)，導(dǎo)致菜心葉部的高積累．

結(jié) 語

本研究以菜心為研究對象，通過對Cd處理前后的葉部、根部的mRNA進(jìn)行高通量測序和分析，獲得了該物種Cd脅迫后的差異表達(dá)基因，這些差異表達(dá)基因的功能和作用機(jī)制還需進(jìn)一步探討．本研究結(jié)果為篩選和培育Cd吸收少、積累低的菜心品種提供了理論基礎(chǔ)．