芳樟葉片轉錄組測序及萜類合成調控相關基因表達分析

曹瑞蘭 胡冬南 周增亮 劉 爽 陳尚钘 劉 娟

（1. 江西農業大學林學院江西省森林培育重點實驗室，江西 南昌 330045；2. 鄭州工業應用技術學院，河南 鄭州 451150）

芳樟（Cinnamomum camphora），又名香樟，是樟科（Lauraceae）樟屬的一種常綠闊葉樹種，主要分布在熱帶和亞熱帶地區[1]。樟樹在建筑、香料、園林造景、生態環境和生態文化建設等方面具有廣泛的用途。根據樟樹精油中主要化學成分不同，將樟樹分為5種類型，包括芳樟（芳樟醇含量58%～92%）、龍腦樟（右旋龍腦含量67%～82%）、腦樟（樟腦含量54%～97%）、油樟（桉葉油素含量32%～52%）和異樟（異-橙花叔醇含量16%～57%）等[2]。其中，從芳樟枝和葉中提取的天然芳樟醇是生產使用最多且需求量最大的一類廣譜性香料，對其特性及其生物合成途徑的研究已成為研究熱點。芳樟醇是芳香類植物香氣的主要揮發性物質，是一種重要的單萜類化合物[3]。芳樟醇合成酶（LIS）是催化單一底物基質牻牛兒基焦磷酸（GPP）轉化成單一產物芳樟醇的關鍵酶。目前，已有研究在不同的芳香植物中開展芳樟醇合成酶基因的表達和克隆等相關研究，但是對芳樟中含量較高芳樟醇合成酶基因差異表達的報道有限。

轉錄組學是從RNA水平研究物種的基因表達，廣泛應用于對特定生物學過程的分子機理研究[4]。因其具有成本低、通量高、時間短、可重復性高等優點，許多植物已經完成了轉錄組測序分析。樟樹轉錄組研究也取得了重要進展。如江香梅等對芳樟醇型、龍腦樟型、樟腦型、油樟型和異樟型等5種樟樹化學類型進行葉轉錄組分析，發現9條Unigene可能參與編碼了芳樟醇合成酶基因[5]。Chen等對芳樟和龍腦樟葉片開展轉錄組測序分析，解析了芳樟醇型和龍腦樟型萜類合成途徑中表達差異的基因，并通過比較分析鑒定了龍腦樟中有關單萜合成酶的3個關鍵基因[6]。已有的文章主要關注樟樹精油中多種萜類化合物合成相關的差異基因，但是他們對芳樟醇含量較高的芳樟中芳樟醇合成酶關鍵基因的挖掘十分有限，需要針對芳樟展開更深入研究。

本研究基于已選育的芳樟優良無性系，并以化學成分較為復雜且芳樟醇含量較低的雜樟為對照組，利用Illumina Hi-seqTM2000高通量測序平臺進行2種化學類型樟樹葉片轉錄組測序。基于獲得的基因序列，進行GO、KOG、KEGG以及Pathway功能注釋，進行差異基因的表達分析，同時重點關注與芳樟醇合成酶合成調控相關的差異表達基因，并開發SSR分子標記。本研究將對推動芳樟分子育種及其樟科植物分子生物學研究和發展具有極大一定的促進作用。

1 材料與方法

1.1 研究材料

本研究材料來源于江西省撫州市金溪縣陸坊鄉石崗村吉源芳樟基地（116°47 E，28°01′ N）的芳樟人工原料林種質資源圃內，選擇長勢良好且無病蟲害的優良芳樟無性系和雜樟實生苗。2017年6月，在優良芳樟無性系3個個體上，分別采集其新葉；在雜樟實生苗的3個個體上采集新葉作為對照。樣品采集后，立即放入液氮中粉碎樣品，-80 ℃低溫保存。芳樟樣品不同個體分別裝樣，實生苗雜樟的葉片則將不同個體混合后勻成3份混樣，進行轉錄組測序。

1.2 方法

1.2.1 RNA提取與文庫構建

采用RNeasy Plant Mini Kit獲得芳樟和雜樟葉（對照）共6個樣品的總RNA，瓊脂糖凝膠電泳檢測其完整性，Nanodrop檢驗總其質量和純度。合格后，磁珠富集和片段化后，用隨機引物（random hexamers）合 成cDNA，AMPure XP beads純化雙鏈cDNA。雙鏈cDNA經末端修復、加測序接頭、加poly（A）等過程，AMPure XP beads篩選片段長度。進行PCR擴增后，經AMPure XP beads純化，最后得到小片段測序文庫，使用Illumina Hi-seqTM 2000進行高通量轉錄組測序及分析。

1.2.2 轉錄組測序與組裝

經CASAVA堿基識別（base calling）分析，使原始測序序列（sequenced reads）被轉化為高通量轉錄測序序列。對所得原始序列（raw reads）進行過濾，去除重復的（N＞5%）、測序質量低的（Q≤20堿基數占50%以上）和帶接頭的reads，得到過濾后序列（clean reads）。采用短reads組裝軟件Trinity[7]將這些短reads連在一起，獲得無法在兩端延伸的序列，稱為Unigene。

1.2.3 Unigene功能注釋

采用BLAST將Unigene與NT（NCBI non-redundant nucleotide sequences）、NR（NCBI non-redundant protein sequences）、Pfam（Protein family）、Swiss-Prot（an annotated protein sequence database）、GO（gene ontology）、COG/KOG（Clusters of orthologous groups/EuKaryotic orthologous groups）、KEGG（Kyoto encyclopedia of genes and genomes）七大數據庫進行對比，獲得Unigene注釋信息。最后以CcEF1α基因作為內參基因，選取潛在的差異基因進行q-PCR，驗證其在不同化學類型中的表達差異。

1.2.5SSRs（simple sequence repeats）檢索與分析

采用MISA對樣品轉錄本數據進行SSR搜索，參數設置分別為：1-10，2-6，3-5，4-5，5-5，6-5（如：1-10表示如果使用單個核苷酸作為重復單元，其重復數≥10才可被檢測到），并且對不同SSR類型在轉錄本的密度分布進行統計。

1.2.4基因表達水平及表達差異量分析

采用RSEM軟件得到比對到轉錄本上的clean reads表達量[8]，對其進行FPKM（Expected number of fragments per kilobase of transcript sequence per millions base pairs sequenced）轉換，并分析基因的表達水平。將所得的轉錄本進行基因差異表達（DEGs）分析。選用的差異基因篩選標準為校準后的P＜0.05且|log2FoldChange|＞1，并對差異表達基因進行GO注釋和KEGG pathway富集分析。

2 結果與分析

2.1 芳樟和雜樟葉精油成分測定

從江西省金溪縣芳樟人工原料林的種質資源圃內采集樟樹葉片樣品，芳樟葉（L-leaf）和雜樟葉（B-leaf）采用GC-MS分析了鮮葉揮發油的化學成分（表1）。在2種化學型的提取物中，芳樟葉的芳樟醇含量達到88.68%，遠高于雜樟葉的芳樟醇含量。

表1 樟樹葉提取物的組成Table 1 Composition of leaf extracts of C. camphora %

2.2 芳香樟和雜樟轉錄組測序與組裝

芳樟葉（L-1 leaf，L-2 leaf，L-3 leaf）和雜樟葉（B-1 leaf，B-2 leaf，B-3 leaf）轉錄組測序后所得raw reads，去除雜質獲得clean reads（表2）。結果表明，質量參數Q20最低為96.17%（芳樟葉），最高達96.91%（雜樟葉）；過濾后不確定的堿基比例N值小于0.01%；其中GC含量為45.72%～46.43%。共生成34.45 GB clean reads數據。

表2 芳樟和雜樟葉測序產量統計Table 2 The statistics sequencing yield of 2 different types of C. camphora

2.3 Unigene注釋與功能分類

利用BLAST獲得的Unigene，分別在Nr、Nt、GO、KOG、Swiss-Prot、Pfam和KEGG數據庫中進行比對，對所注釋的Unigene進行統計，并通過BLASTN將其比對到Nt（E值＜0.000 01），獲得其功能注釋信息和具有最高相似性的蛋白（表3）。結果顯示，312 457條Unigene在NR數據庫中存在98 579條同源匹配到信息，Nt 67 148條，Ko 40 723條，SwissProt 80 310條，Pfam 85 151條，GO 87 053條，KOG 33 177條。

表3 Unigenes注釋到不同數據庫統計Table 3 The statistics of Unigenes annotation to different databases

2.4 Unigene的GO功能注釋和分類

基于Nr的注釋信息，在GO數據庫對芳樟和雜樟的葉轉錄組Unigene進行生物學特征功能分類。結果顯示，GO數據庫注釋Unigene可分為生物過程（biological process）、分子功能（molecular function）和 細 胞 組 分（cellular component）3個GO類別，56個小組。在生物過程中，共含有207 529個GO條目，25個功能組，其中“代謝過程”45 038個，“細胞進程”47 138個，“單一生物體”過程35 746 個。在細胞組分中，共涉及132 394個GO條目，21個功能組，其中“細胞”26 278個，“細胞要素”26 259個，“細胞器”17 788個。在分子功能中，共涉及95 197個GO條目，分為10個功能組，其中“催化活性”有36 459個，“結合活性”41 049個，“細胞聚集”2個，“細胞外基質組分”7個。

2.5 差異表達基因的表達分析

為了分析芳樟（L-leaf）和雜樟（B-leaf）2種不同化學類型葉片的差異基因表達，參考edgeR進行差異基因分析，篩選閾值選擇為false discovery rate（FDR）≤0.05，|log2FoldChange|≥1，共鑒定出11 685個單基因，其中芳樟相對于雜樟有6 481個上調基因和5 204個下調基因，芳樟特有差異基因794個，雜樟特有842個差異基因，10 049個基因在2種化學類型中均有表達（圖1）。

圖1 樣品間差異表達基因的數量Fig. 1 Number of differentially expressed genes between 2 samples

2.5.1 差異基因GO功能注釋

將芳樟和雜樟葉所獲得的轉錄本用DEGseq軟件進行比較分析，獲得差異表達基因11 685個，將差異基因進行GO分類統計（表4），參與生物學過程差異基因分為20類，其中磷代謝過程（727個）、含磷化合物代謝過程（725個）和大分子修飾（640個）中所占比例最高。參與分子功能中包含了20個功能小類，差異基因主要分布在催化活性（2 810個）、轉移酶的活性（1 180個）和核苷酸綁定（821個）中。在細胞組分這一類中，外部封裝結構（98個）和細胞壁（79個）中差異基因較高。

表4 差異基因在GO中的功能注釋Table 4 Functional annotations of different genes in GO

2.5.2 差異基因KEGG富集分析

將獲得的11 685個DEG映射到KEGG代謝途徑數據庫，共有2 200條DEG獲得功能注釋，占DEG總數的18.83%，其中歸入KEGG的代謝途徑有116條。將富集的上調和下調差異基因最多的前20條KEGG代謝途徑列出（表5），主要包含植物－病原物互相作、糖代謝、植物激素信號轉導、內質網蛋白加工、苯丙烷生物合成和吞噬體等，其中，上調差異基因中涉及倍半萜和三萜生物合成的Unigene有12個，下調差異基因中發現單萜類生物合成的Unigene有5個。

2.5.3 與萜類生物合成調控相關差異表達基因

為了探索不同萜類化合物積累模式的調控機制，對萜類生物合成相關基因的表達進行深入分析。共鑒定39個編碼萜類合成酶的Unigenes，其中單萜生物合成有6個，二萜生物合成有7個，倍半萜和三萜生物合成有12個，萜類骨架合成14個。表5是萜類生物合成調控的差異基因，在雜樟葉中有5個單萜合成酶（TPS14）基因、3個二萜合成酶基因、2個萜類骨架合成酶基因顯示上調，而在芳樟葉中有1個單帖生物合成基因、4個二萜合成酶基因12個Unigene參與倍半萜和三萜類合成酶基因和12個萜類骨架合成酶在表達上調（表6）。

表5 前20個上調和下調的差異基因KEGG代謝通路顯著性富集列表Table 5 The top 20 up-regulated and down-regulated differential genes KEGG metabolic pathways with significant sexual enrichment

表6 與萜類生物合成調控相關差異表達基因Table 6 The list of differentially expressed genes related to the regulation of terpenoid biosynthesis

2.5.4 芳樟醇合酶（TPS14）基因的表達量分析

為驗證RNA-Seq分析獲得的萜類生物合成基因的表達模式，本研究采用q-PCR方法檢測了芳樟（L-1leaf、L-2leaf和L-3leaf）和雜樟（B-1leaf、B-2leaf和B-3leaf）2種化學類型中3個單帖合成酶（TPS14）基因的表達水平（圖2a）。結果顯示，3個TPS14基因在2個化學類型中的q-PCR表達模式與RNA-Seq（圖2b）的FPKM的結果基本一致，均表現為3個TPS14基因在芳樟中的表達量比雜樟中的較低。

圖2 TPS14基因在兩種化學型中葉表達模式分析Fig. 2 Analysis of the expression patterns of TPS14 genes in 2 chemical types of leaves

2.6 SSR分析情況

利用MISA對轉錄組測序獲得的312 457條Unigene進行檢索，共找到符合條件的SSR位點有134 851個，分布于88 016條Unigene中。其中含單個SSR位點的Unigene有103 478條，含2個及2個以上SSR位點的Unigene有31 373條。復合型SSR數目為20 277個（表7）。

表7 轉錄組中SSR信息分析Table 7 SSR information analysis in transcriptome

3 結論與討論

芳樟的枝葉經減壓精餾分離而來的高純度芳樟醇，在香精、香料和醫療保健等方面具有廣泛的應用價值[9]。而雜樟，通常是指野生樟樹中精油主成分不明顯的一類，因在自然演化歷史中，經反復雜交和環境因素的影響，其化學成分和含量十分復雜，生產利用途徑有限[10]。目前關于樟樹的研究主要集中在生長發育[11-12]、栽培管理[13-14]、油脂組成及提取工藝[15-16]等方面，有關其分子生物方面研究有限，尤其是其主要成分合成途徑及其分子機理尚不清楚。本研究通過Illumina HiSeqTM 2000高通量基因測序技術對芳樟和雜樟的葉轉錄組進行比較分析，轉錄組測序共生成34.45 GB clean數據，組裝成312 457個Unigene，鑒定出11 685個差異表達Unigenes。相對雜樟，芳樟葉中有6 481個上調的Unigenes和5 204個下調Unigenes，其中有2 200條差異基因在KEGG代謝途徑數據庫中獲得功能注釋，占DEG總數的18.83%，被歸入116條KEGG的代謝途徑中。

KEGG pathways分析結果表明，注釋到次生代謝生物合成途徑中的Unigene有1 184條，占2.91%。樟科的各個家族中萜類化合物含量豐富[17-18]，目前已在山蒼子（Litsea cubeba）[19]、細毛樟（Cinnamomum tenuipilum）[20]、土肉桂（Cinnamomum osmophloeum）[21]、月桂（Laurus nobilis）[22]等類群中成功鑒定參與萜類生物合成的關鍵基因。針對樟樹的不同化學類型，研究者也展開了萜類合成酶關鍵基因的探索。如江香梅等通過樟樹5種化學型的轉錄組分析，鑒定發現了參與編碼芳樟醇的9個關鍵基因，且這些基因在芳樟中表達水平遠高于油樟等其他化學類型[5]。這與本研究的結果相似，本研究發現的39個與萜類生物合成調控相關基因，涉及單帖合成途徑的的Unigene有6個，二萜合成通徑的Unigene有7個，倍半萜合成通徑的Unigene有12個，萜類骨架合成通徑的Unigene有14個。其中，5個單萜合成酶（TPS14）基因、3個二萜合成酶（E1.14.11.13）基因、2個萜類骨架合成酶（BCHP和FOLK）基因在芳樟葉中表達水平較低，而在雜樟葉中表達水平較高。

江香梅等對芳樟和龍腦樟的轉錄組分析中發現它們存在67條與萜類合成相關的差異表達基因，且多個基因被注釋為同一種酶[5]。本研究基于芳樟和雜樟葉轉錄組數據，鑒定了39個萜類合成相關的差異表達基因。造成結果的差異，可能是由于與所選用材料的萜類物質含量不同。同時也發現5個TPS14基因均被注釋為芳樟醇合酶。Chen等的結果顯示，TPS14在芳樟中表達量較低，且在樟樹的5種化學型葉片熒光定量表達分析發現相同的規律，即TPS14在芳樟醇含量較高的芳樟中反而表達最低[6]。本研究的TPS14基因轉錄組FPKM數據和q-PCR表達量也證實了該結論，潛在的3個TPS14基因在雜樟的表達量要高于芳樟。然而，在擬南芥中，TPS14被認為是花部器官中調控芳樟醇合成的關鍵酶[23]。已有報道發現，TPS少數氨基酸的變化就會導致其特定酶合成萜類譜發生劇烈變化[24]。Liu等通過研究油茶的萜類合成代謝途徑，發現選擇性剪切可改變關鍵酶的轉錄，進而影響芳樟醇合成[25]。TPS14基因在芳樟中的轉錄過程及其功能需要進一步探索和證實。

目前，在樟樹精油產業中，種質資源混亂和缺少優良無性系是制約其發展的主要因素之一。由于不同化學類型的精油含量及其品質很難依靠樟樹外部形態進行精準劃分，而精油提取和化學成分分析的方式存在時效慢、受季節和生長發育階段的影響等缺點，因此急需新的鑒定手段。SSR分子標記法，普遍認為具有共顯性、重復性好、多態性豐富、操作簡便等優點[26]，目前已成為開展植物品種鑒定、遺傳多樣性分析、遺傳圖譜構建和分子輔助育種[27]的主要方式之一。張震等以油茶轉錄組測序數據為基礎，利用MISA軟件對SSR進行搜索，發現104 515個SSR，分布在80 724條Unigene中[28]。本研究基于芳樟和雜樟葉的轉錄組數據SSR位點分析，有88 016條Unigene鑒定到SSR位點。SSR特征分析的結果，為芳樟和雜樟及其同屬物種的目標基因的標定、遺傳多樣性分析、遺傳圖譜構建和基因組差異分析提供了數據依據。