鎘脅迫后旱柳轉錄組變化分析

曹繼敏 李雙財 何德

（西南林業大學生命科學學院，昆明 650224）

隨著工業的發展，鎘（Cd）成為全球污染最主要的重金屬之一，如工業廢氣會將Cd帶到氣流層中，并以大氣沉降、降雨以及下雪等方式進入土壤中，造成土壤Cd污染；同時，農業上的污水灌溉以及含Cd藥劑的使用都會使Cd進入土壤，從而加重土壤的Cd污染。此外土壤中的Cd主要通過根的吸收進入植物體內，然后與根細胞的細胞壁結合，后經木質部運輸到莖、葉、果實等器官中積累，抑制植物的生長發育［1］。Cd可以通過食物鏈進入動物和人體內，對機體造成嚴重傷害。過量Cd積累導致肺紋理增多、紊亂而模糊，白細胞和中性粒細胞數量增多，引起肺泡Ⅱ型上皮細胞毛細血管內皮細胞受損，使肺泡間隔增厚，間質纖維增生，也會直接作用于骨骼，使有機體骨骼礦化發生障礙、骨鈣溶出的增加以及骨膠原和骨的固化作用異常等［2-3］。20世紀60年代發生在日本神通川流域的“骨痛病”，原因就是當地居民食用富含重金屬Cd的大米造成的。在我國也發生了由于土壤Cd污染日益嚴重出現了例如“鎘大米”，農作物減產等一系列惡性事件［4］。

考慮到傳統的物理、化學修復重金屬污染的方法中存在不可避免的缺陷［5］，植物修復以其成本低、安全可靠、對環境干擾小且改善生態環境等優點而備受關注［6］，木本植物因為其具有龐大的根系、巨大的生物量、發達的維管組織、堅固的木質組織、以及高效的蒸騰作用，并且可以有吸收多種重金屬的特點而成為環境重金屬污染治理的重要方向［7］。旱柳（Saliz matsudana）作為一種常見的木本植物，具有生長快、易繁殖、生物量大、根系發達、再生能力強和對多種重金屬耐受的特點［8］，此外旱柳的葉與根等新陳代謝較快的器官可以大量積累重金屬Cd，而其它的主要營養儲存器官（果實、籽粒、塊莖等）對于重金屬積累較少，因此旱柳具有巨大的治理重金屬污染的潛力［9］。

當前對旱柳重金屬耐受性的研究主要集中在對重金屬的植物體內富集情況的調查、影響因素、生理變化等方面，對重金屬富集后植物耐受的分子機理研究較少，尤其在轉錄水平上研究很少。Claudia等［10］對柳樹耐Cd無性系葉片中Cd的沉積研究發現Cd主要分布于嫩葉的葉尖以及老葉的葉頸，水培時間較短的旱柳主要以嫩葉為主。目前我國土壤一般Cd污染濃度為2.5 mg/L，徐愛春［11］研究表明Cd濃度為50 mg/L時對旱柳具有較大影響效應，本研究以旱柳葉片為實驗材料從轉錄組水平著手，探究旱柳在Cd脅迫條件下的耐受性反應情況，為旱柳修復土壤Cd污染提供理論指導。

1 材料與方法

1.1 材料

本研究所用的旱柳枝條取樣自昆明市北京路延長線的一株旱柳。一年生旱柳枝條被截取長度為25 cm繼而用改良后霍格蘭德營養液［12］進行水培。水培后的旱柳隨機分成3組（每組5盆，每盆3棵），1組為對照組，其余2組為試驗組。Cd（NO3）2作為唯一的Cd來源對試驗組進行脅迫處理，使營養液中Cd濃度分別為2.5 mg/L和 50 mg/L，營養液每3 d換一次。研究材料為3組的第0、1、7與30天并處于同一位置的葉片，采樣方法使用隨機混采法；所有試驗葉片采集完成后液氮急速冷凍后放置于-80℃冰箱內保存備用。對照組用CK表示，2.5 mg/L處理組用L表示，50 mg/L處理組用H表示；樣品編號如下:第0天為 CK0，第1 天分別為CK1、L1、H1，第7天分別為CK7、L7、H7，第30天分別為CK30、L30、H30。

1.2 方法

1.2.1 總RNA提取及樣品檢測 總RNA提取使用OMEGA公司的E.Z.N.ATM Plant RNA Kit試劑盒進行提取，提取方法嚴格按照試劑盒說明書進行。并將獲得的總RNA用瓊脂糖凝膠電泳分析RNA的降解程度及是否存在污染，用液相色譜Agilent 2100和OD260/280的比值檢測RNA的純度。旱柳總RNA的提取、質量檢測及無參考基因組序列的轉錄組高通量測序分析均由深圳華大基因生物科技有限公司進行。

1.2.2 文庫構建及轉錄組測序 樣品檢測合格后，加入Fragmentation buffer將mRNA打斷成短片段，用六堿基隨機引物（Random hexamers）以mRNA為模板合成一鏈cDNA，然后加入緩沖液、dNTPs和DNA polymerase I合成二鏈cDNA，雙鏈cDNA再利用AMPure XP beads純化。純化后雙鏈cDNA末端修復、加A尾并連接測序接頭，片段大小選擇是用AMPure XP bead，最終的cDNA文庫通過PCR富集得到。cDNA 文庫質量用Agilent 2100檢測，再通過Illumina Hi-Seq2500平臺完成RNA-Seq測序分析。

1.2.3De nove組裝 獲得的原始數據Raw reads經過濾得到Clean reads。使用Trinity對clean reads組裝得到contigs，之后再用Tgicl對轉錄組樣本進行兩次聚類去冗余得到最終Unigene。Unigene分為兩部分，一部分是clusters，是進一步去冗余后的結果；其余的是singletons，指沒有聚類、單獨的 Unigene。每個樣品測序產出不少于 6 Gb Clean data，Q30堿基百分比達到85%。

1.2.4 生物信息學分析 使用BLAST 對Unigene序 列 與 Nr、Nt、KOG、KEGG、GO、Swiss-Prot及InterPro數據庫進行比對（Evalue＜1e-5），獲取與旱柳Unigene 具有最高序列相似性的蛋白，得到Unigene 的蛋白功能注釋信息。

本研究從兩個角度分析旱柳在Cd脅迫后的差異基因表達情況:（1）同一時間不同濃度Cd脅迫的差異基因表達情況的對比分析，命名為L1-CK1（即L1和CK1兩個樣品的基因表達對比，下同）、H1-CK1、L7-H7、L30-CK30、H30-CK30；（2）同一濃度不同時間的差異基因表達情況對比分析，分為2.5 mg/L組和50 mg/L組，命名為L1-CK0、H1-CK0、L7-L1、H7-H1、L30-L7、H30-H7。

根據差異表達分析結果，使用 Blast2GO和WEGO軟件對差異表達的Unigene進行GO功能分類［13］。將差異表達的Unigene注釋到KEGG數據庫中得到注釋信息，之后進行pathway分析［14］。得到2.5 mg/L和50 mg/L濃度Cd脅迫后旱柳的生理生化響應結果。

1.2.5 差異表達基因q-PCR驗證 利用 primer premier 5.0 對Cd脅迫響應的相關基因和內參基因肌動蛋白基因（ACT）進行 q-PCR 引物設計。在BIORAD公司的C1000TM熒光定量PCR儀進行real-time PCR反應，每個樣品重復3次，反應體系20 μL。以內參基因ACT為對照，利用公式2-△△Ct計算其相對表達量。

2 結果

2.1 總RNA質檢結果

10個樣品中除CK7外9個樣品的RNA質量較好，能滿足建庫測序要求，可用于RNA轉錄組測序，且總量滿足1次但不足2次建庫需要。CK7因RNA提取質量不佳，所以沒有進行轉錄組測序。

2.2 轉錄組測序與De nove組裝

9個樣品的raw reads在72.22 Mb到85.59 Mb之間，過濾后總clean reads在 65.17 Mb到67.15 Mb之間；Q20、Q30以及過濾前后reads的比率分別是95.04%到 97.38%，86.40%到 90.91%，72.75%到89.4%（表1），這表明低質量的堿基比率低，測序質量好，可以用于de novo組裝。用Trinity對其進行組裝，所有樣品的總reads均超過10萬，總長度超過7千萬，平均長度超過600，N50超過1 000，GC含量42%左右。之后使用 Tgicl 對轉錄本進行聚類去冗余得到Unigene，聚類后的 Unigene質量（表2），10個樣品的Unigene數為53 705到62 933個，總長度從43 247 382 nt 到50 853 969 nt，N50 值都在1 200之上，這表明組裝質量高，可用于后續分析。

表1 過濾前后的reads質量統計表

表2 Unigene的質量指標

2.3 對Unigene的功能注釋

得到的Unigene有83.93%在7個數據庫（Nr、Nt、KOG、KEGG、GO、Swiss-Prot及 InterPro）中任意一個注釋到。Nr數據庫進行比對篩選，篩選條件為E值小于10-5，結果顯示:有74 919條Unigenes在Nr數據庫中找到同源蛋白序列，占總Unigenes數的73.02%。Nt數據庫注釋到的Unigenes最多，占總Unigenes的81.12%，但被全部數據庫注釋到的Unigenes僅為29.04%，說明有許多旱柳Unigenes的功能還沒有完全明確。

2.4 Unigene GO功能分類

根據 Nr 數據庫注釋的信息結果顯示，有51 494條Unigene 映射到GO不同的功節點（Term）上，根據生物過程（Biological process）、細胞組成（Cellular component）、分子功能（Molecular function）進行分類。由于經常有同一個轉錄本映射到不同節點現象，所以在生物過程中有106 517條，占42.69%，參與代 謝 過 程（Metabolic process）（22 955，21.54%）、細 胞 過 程（Cellular process）（21 034，19.75%）的Unigene最多；在細胞組成中有92 900條，占37.24%， 參 與 膜（Membrane）（15 182，16.34%）、細胞器（Organelle）（13 010，14.01%）、細胞（Cell）（16 820，18.11%）和細胞部分（Cell part）（16 671，17.95%）的Unigene最多；在分子功能中有50 054條，占20.06%，參與催化活性（Catalytic activity）（20 318，40.59%）和結合蛋白（Binding）（20 533，41.02%）的Unigene最多（圖1）。

圖1 Unigene的GO分類圖

2.5 Unigene 的KEGG注釋和代謝通路分析

將旱柳差異表達基因在KEGG數據庫中進行比對分析顯示，53 966個Unigene分別注釋到KEGG數據庫的6個一級層級和21個二級層級中，共包含137條通路。其中在同一個pathway通路中10個以上差異表達基因主要涉及細胞過程、環境信息處理、遺傳信息處理、人類疾病、代謝以及有機系統6個方面。這6個方面中代謝類的涉及最多，包括碳水化合物、氨基酸、多糖生物合成、能量、次生物質等；涉及人類疾病的最少，只有2個，分別是內毒素和代謝疾病（Endocine and metabolic diseases）和耐藥性∶抗菌素（Drug resistance∶Antimicrobial）。在旱柳Cd脅迫后KEGG分析中發現涉及人類疾病的基因變化，這也表明Cd對旱柳中有關人疾病的內毒素和抗菌性的代謝有影響。同時，pathway富集結果也表明，前述的25個Cd脅迫響應候選基因主要參與代謝途徑（Metabolic pathway）（10541，19.53%）、次生代謝產物的合成（Synthesis of secondary metabolites）（5473，10.14%）、植物MAPK信號通路（Plant MAPK signaling pathway）（1592，2.95%）、 植 物 病原互作（Plant pathogen interaction）（1761，3.26%）、植物激素信號轉導等通路（Plant hormone signal transduction）（1862，3.45%）（表3）。由此推測旱柳可能主要是通過這些代謝通路調控Cd抗逆性過程。

表3 富集KEGG通路分析

2.6 差異表達基因分析

FDR（False discovery rate）值越小，表明基因表達差異越大，差異越顯著。本研究在差異基因時設置的閾值為Fold Change≥2.00，FDR≤0.001，Q-value＜0.05。

由表4可知，相同時間不同濃度的差異表達基因（Differentially expressed genes，DEG）分析顯示:第1天組（L1-CK1和H1-CK1）中的差異基因數為8 550個（3 660個上調和4 890個下調）；第7天組（L7-H7）的差異基因數為3 444個（1 307個上調和2 137個下調）；第30 天組（L30-CK30和H30-CK30）中的差異基因數為11 428個（5 055個上調和6 373個下調）。相同濃度不同時間的DEG分析顯示:2.5 mg/L組（L1-CK0，L7-L1，L30-L7）中的差異基因數為26 623個（14 531個上調和12 092個下調；50 mg/L組（H1-CK0，H7-H1，H30-H7）中的總差異基因數為32 154個（18 318個上調和13 836個下調），這說明50 mg/L濃度的Cd脅迫比2.5 mg/L的對旱柳的影響更加明顯。此外，在L1-CK1、H1-CK1、L30-CK30和H30-CK30這4組中有共同的170個DEG（其中46個共同上調和29個共同下調），它們的基因表達量差異達2倍（|Log2ration|≥1）以上。170個DEG中又有25個其表達量差異是大于1 000倍（|Log2ration|≥10）以上的，其中19個是基因表達上調的，6個是基因表達下調的（表5）。這25個基因表達變化如此之大，說明它們在旱柳受Cd脅迫耐受過程中等起到了重要作用，可以作為旱柳Cd脅迫響應的重要候選基因。

表4 差異表達基因統計

表5 抗逆相關基因詳細情況

2.7 差異基因的q-PCR

黃酮類化合物生物合成和油菜素內酯合成的關鍵基因:FLS，F3H和ROT3，這3個基因在旱柳Cd脅迫響應代謝通路中表達量顯著上升，鈣依賴性蛋白激酶2基因表達量顯著下調，因此我們選擇這4個基因進行q-PCR驗證轉錄組測序分析的準確性和可靠性，肌動蛋白基因（ACT）因其具有穩定性而被選擇為內參基因［15］。圖2-6的q-PCR結果顯示:5個差異表達基因在轉錄組測序中的表達變化與熒光定量結果基本一致，說明轉錄組測序分析的結果真實可靠。

圖2 CL11387.Contig1基因q-PCR結果

圖3 Unigene3229_All基因q-PCR結果

圖4 CL13372.Contig3基因q-PCR結果

圖5 CL1644.Contig4基因q-PCR結果

3 討論

圖6 CL10874.Contig2_All基因q-PCR結果

旱柳是中國土生土長的高大喬木，具有易扦插、易成活、生長快、枝葉發達和抗逆性強等特點，是一種非常適合用修復土壤重金屬污染的物種。為研究旱柳響應不同濃度Cd脅迫的基因表達情況，選擇Cd脅迫濃度分別為2.5 mg/L和50 mg/L的旱柳葉片作為轉錄組測序樣品，通過旱柳基因表達的差異來研究旱柳響應Cd脅迫的分子機制。轉錄組數據進行過濾、拼接、組裝后共獲得102 595個基因。

近年來對植物響應重金屬脅迫后的基因表達研究已取得一系列成果，鑒定和推測出不少關于重金屬響應的基因，包括ABC轉運家族、鋅和錳轉運蛋白家族、重金屬ATP酶（HMAs）家族、陽離子擴散促進者（CDF）家族、植物螯合肽合成基因（PCS）、金屬硫蛋白基因（MT）和金屬忍耐蛋白基因（MTP）等［16］。本研究發現在旱柳響應Cd脅迫的過程中，金屬硫蛋白基因、ABC轉運家族、鋅和錳轉運蛋白家族的基因表達會隨著Cd濃度的增高和脅迫時間的增加而有較大的變化，L7-L1、L30-L7、H7-H1、H30-H7組中金屬硫蛋白共9個基因表達上調、11個下調；ABC轉運家族共90個上調、88個下調，其中H30-H7中32個上調，50個下調；鋅和錳轉運蛋白家族有28個下調，8個上調，上面列的其它重金屬響應基因的差異表達在本研究沒有發現。這些結果表明旱柳在應對Cd脅迫的過程中金屬硫蛋白基因、ABC轉運家族及鋅和錳轉運蛋白起了很大的作用，也說明Cd脅迫會影響旱柳對其他金屬離子的吸收，鋅和錳轉運蛋白表達量的變化可能是因為Cd在生命體內會置換鋅，導致生命體需鋅的蛋白質發生“饑餓”，為了生命正常活動的進行，旱柳加快了鋅的轉運過程，從而使鋅、錳蛋白轉運蛋白家族的表達發生變化。ABC轉運家族基因表達量發生明顯變化，說明旱柳在Cd脅迫響應中細胞物質跨膜轉運主要是通過ABC轉運家族蛋白來實現的，也說明ABC轉運家族在旱柳細胞運送重金屬相關物質中起重要作用。

黃酮類化合物是植物體內一類在抗菌、抗逆性上具有重要的作用的化合物，是一種良好的抗氧化劑，能夠清除生物膜周圍和細胞內的H2O2［17-18］。有

研究表明，黃酮化合物的代謝通路在水稻響應Cu2+和Cd2+的脅迫以及組培旱柳苗響應Cd2+中起了積極作用［19-20］。孫愛清等［21］利用 Solexa高通量測序對花生干旱響應基因表達基因分析，發現有9個類黃酮代謝相關基因在干旱脅迫下顯著表達；宋中邦［22］通過煙草基因組數據庫中搜索獲得兩個類黃酮合成酶基因，即NtFLS1和NtFLS2，二者表達模式一致，均在葉片中高水平表達。本研究通過KEGG代謝通路圖進行分析發現，旱柳黃酮類化合物的生物合成通路在樣品中均有大量的差異基因上調，L1-CK1、H1-CK1、L30-CK30、H30-CK30組 中 一 共 170個差異基因，其中有46個共同上調基因，29個共同下調基因，而這46個共同上調基因中包含了FLS、F3H基因。這說明旱柳在應對Cd脅迫時，黃酮類化合物起了重要作用。

同時，本研究發現油菜素內酯合成通路的 3，6-脫氧油菜素淄酮酶（ROT3）在除L1-CK1外的其他樣品中都是上調。肖瑞雪等［23］研究發現油菜素內酯是一種被譽為第六激素的天然激素，其生理活性甚至超過現在5種激素，能提高植物的抗逆性。靳開川等［24］結合遺傳學、基因與蛋白質組學、細胞生物學等多學科方法和手段，發現油菜素內酯在植物的抗逆性（抗干旱、高鹽、高溫、低溫、重金屬）過程中的起重要作用。王喆等［25］研究表明油菜素內酯可能會與膜蛋白結合，然后通過減輕重金屬脅迫來提高代謝活性，從而降低植物對重金屬的攝取能力［26］。本研究在 H1-CK1、L30-CK30、H30-CK30中油菜素內酯合成通路的 3，6-脫氧油菜素淄酮酶（ROT3）上調是由于Cd2+在旱柳葉片中積累，旱柳為減輕Cd2+帶來的損害，而促進了油菜素內酯的合成。

本研究對旱柳受Cd脅迫的轉錄組進行了差異基因的分析，找到了受脅迫的相關基因以及與之相關的代謝通路，但這些基因是如何互相作用調節旱柳Cd脅迫，還需要在分子方面進行更加深入的分析。

4 結論

旱柳在受重金屬Cd脅迫后展現出一定程度的脅迫抗性，同時基因表達與代謝過程發生了重大的變化。（1）鎘脅迫后不同濃度不同時間的旱柳轉錄組測序共獲得102 595個Unigenes；（2）篩選得到與Cd脅迫響應密切相關的基因25個，其中金屬硫蛋白、ABC轉運蛋白、鋅和錳轉運蛋白基因在旱柳植物抗逆性過程中發揮重要作用；（3）Cd脅迫的GO條目主要集中在代謝過程、細胞過程、膜、細胞器、細胞、細胞部分、催化活化和結合蛋白上；（4）油菜素內酯合成通路的3，6-脫氧油菜素淄酮酶（ROT3）和黃酮類化合物合成通路的黃酮醇合成酶（FLS）、黃烷酮-3-羥化酶（F3H）均明顯上調。