王芬 裴會敏 文狄 陳志 劉榮 姚玉仙 馬媛

摘 要：為探究茶樹中茶多酚等產物代謝途徑的相關基因，該研究以貴州都勻地區(qū)福鼎大白種茶樹的根莖葉為對象，利用高通量測序技術構建茶的轉錄組數(shù)據(jù)庫并篩選其根莖葉差異表達基因。結果表明：共獲得70.88 Gb Clean Data，各樣品Clean Data均達到6.33 Gb，Q30堿基百分比在93.22%以上。將Clean Reads與中國種茶樹參考基因組進行序列比對，比對效率從87.83%到91.14%。基于比對結果，進行可變剪接預測分析和基因結構優(yōu)化分析，發(fā)掘新基因13 531個，其中10 244個得到功能注釋。利用FPKM進行基因表達量分析，根據(jù)基因在不同樣品中表達量識別差異表達基因。葉與莖的差異基因有5 595個，其中2 769個在莖中上調，2 826個下調，葉與根有9 650個差異基因，5 056個上調，4 594個下調，莖與根中有5 644個差異基因，2 938個上調，2 706個下調，并通過GO和KEGG分析，將差異基因進行功能注釋和富集分析。上述結果為揭示都勻地區(qū)福鼎大白種茶參與類黃酮、茶氨酸和咖啡堿等代謝途徑相關的基因提供了參考，為選育優(yōu)良品種等提供了理論依據(jù)。

關鍵詞：都勻地區(qū)，福鼎大白種，根莖葉，高通量測序，差異基因

中圖分類號：Q945

文獻標識碼：A

文章編號：1000-3142（2020）09-1269-12

Abstract：Tea tree is rich in catechins，theanine，caffeine and other metobolite of health fuction. In order to study the related genes of the metabolisms of the polyphenols. We use high-throughput sequencing technology to study the root，stem and leaf of Fudingdabai tea and find differential expression genes （DEGs）. The results showed that 70.88 Gb Clean Data was obtained，6.33 Gb Clean Data is in each sample and Q30 is more than 93.22%. We map the Clean Reads to reference genome，the blast result is from 87.83% to 91.14%. Then，alternative splicing and gene structure optimization was analyzed. There are 13 531 new genes，in which，10 244 genes were annotated. GO and KEGG functional annotation and enrichment analysis were carried out in differential expression genes，which were identified according to gene expression level in different samples. There were 5 595 DEGs between leaf and stem，2 769 genes were up-regulated and 2 826 genes were down-regulated. 9 650 DEGs were found beween leaf and root，5 056 genes were up-regulated and 4 594 genes were down-regulated. 5 644 DEGs between stem and root，2 938 genes were up-regulated and 2 706 genes were down-regulated. The results are expected to provide reference for recognizing genes of catechins，theanine，caffeine pathways，provide the theoretical basis for breeding improved seeds.

Key words：Duyun，F(xiàn)udingdabai，root，stem，leaf，high-throughput sequencing，differential expression gene

都勻地處云貴高原苗嶺山脈南側，具有低緯度、高海拔、寡日照、多云霧的自然生態(tài)環(huán)境，適宜茶樹生長，資源豐富（陳世軍等，2017）。不同品種的適制性、抗逆性和生長勢不同，從而導致茶品質和產量的不同。20世紀90年代以來，都勻市等地開始大量引種單產高、品質好、適應性強的福鼎大白種茶樹，增加了茶樹的遺傳多樣性。目前，關于都勻地區(qū)的茶種質資源保護和抗性研究薄弱（張麗娟等，2016），特別是在轉錄組層面上的研究較少，在一定程度上限制了茶產業(yè)的應用發(fā)展。福鼎大白種茶園占全市茶園面積85%以上，因此開展都勻地區(qū)福鼎大白種根、莖、葉的轉錄組研究對茶樹在品質、質量及抗逆相關的遺傳改良上具有重要意義。

由于高通量測序技術能快速全面獲取研究對象在某一狀態(tài)下基因轉錄信息，因此該技術廣泛應用于生物體轉錄組基因表達分析，能準確發(fā)掘重要功能基因。Shi et al.（2011）用RNA-seq技術對龍井43的嫩葉、成熟葉、莖、幼根、花蕾以及成熟種子進行測序，總共獲得127 094個基因，其中包括與茶特有香氣相關的代謝產物如咖啡堿、茶氨酸、黃酮等合成相關的基因。王君雅等（2019）利用轉錄組技術對龍井43和中茶126進行差異基因表達分析，結果表明99個表達模式完全相反的基因參與的生物過程主要有MAPK信號通路、谷胱甘肽和苯丙烷等。利用轉錄組測序技術可識別編碼次生代謝的調節(jié)基因和結構基因，也可用來預測未知基因的功能，是進行基因組功能研究的重要工具（Novaes et al.，2008）。本研究利用高通量測序技術，對都勻地區(qū)福鼎大白種茶樹根莖葉中類黃酮相關代謝途徑進行研究，探索根莖葉中差異表達基因參與的代謝途徑，為探討都勻地區(qū)福鼎大白種茶樹生長發(fā)育及組織間功能差異的分子機制提供理論依據(jù)。

1 材料與方法

1.1 材料

供試材料選用黔南州都勻地區(qū)種植的茶樹-福鼎大白種。茶苗選用黔南民族師范學院茶園的扦插苗。選取9株培養(yǎng)條件相同長勢相同的茶苗分為3組，每組3株，在第1組中取嫩根、嫩莖和嫩葉分別放入3個離心管作為第1組生物學重復的樣本，分別標為根的第1個生物學重復，莖的第1個生物學重復，葉的第1個生物學重復。每個樣本3個生物學重復，剩余的兩組生物學重復分別在剩下的兩組茶苗中取樣，共9個樣本，分別為葉1、葉2、葉3、莖1、莖2、莖3、根1、根2和根3。采樣時在茶樹上剪取發(fā)育階段相似、部位相同、充分伸展的葉片、莖和根，樣品隨后迅速放入液氮罐中速凍后放入-80 ℃冰箱中儲存?zhèn)溆谩?/p>

1.2 轉錄組測序

首先將采集的福鼎大白種茶苗根莖葉樣品放在干冰中送往北京百邁客生物科技有限公司進行二代Illumina高通量測序，對測序數(shù)據(jù)進行質量控制，并將二代得到的轉錄本與已測序的中國種茶樹基因組通過序列比對的方法進行比較。然后進行轉錄組文庫質量評估、SNP/InDel分析、可變剪接事件預測、差異基因分析和差異外顯子分析等。

Illumina二代測序數(shù)據(jù)提交至NCBI的SRA數(shù)據(jù)庫，BioProject的編號為PRJNA562747。

2 結果與分析

2.1 測序數(shù)據(jù)質量分析

從茶苗上分別取福鼎大白的根莖葉作為實驗材料進行轉錄組測序，堿基質量值達到Q30的在93.22%以上，GC含量為44.56%～45.85%（表1），結果表明轉錄組文庫質量符合分析要求。

2.2 轉錄組與中國種茶樹基因組比對

2.2.1 比對效率 本研究利用Hisat2（Kim et al.，2015）對測序數(shù)據(jù)與中國種茶樹基因組（Wei et al.，2018）進行比對，利用StringTie（Pertea et al.，2015）對比對上的Reads進行組裝和定量，比對到參考基因組上的Reads占Clean Reads的百分比為87.83%～91.14%（表2）。

2.2.2 SNP/InDel分析 利用GATK（Mckenna et al.，2010）軟件對Hisat2比對結果中的SNP位點和InDel進行識別，進而分析基因表達水平和基因功能，SNP分析表明基因區(qū)SNP位點數(shù)多于基因間區(qū)，轉換型SNP多于顛換型SNP（表3）。InDel主要存在于內含子區(qū)和基因間區(qū)（圖1）。

2.2.3 可變剪接預測 基因通過轉錄生成前體mRNA，再經過不同的剪接，產生不同的成熟mRNA，翻譯為不同的蛋白質。利用Asprofil（Florea1. 基因間; 2. 基因內; 3. 內含子; 4. 上游; 5. 下游; 6. 受體剪切位點; 7. 供體剪切位點; 8. 剪切位點; 9. 起點缺失; 10. 移碼; 11. 密碼子缺失; 12. 密碼子插入; 13. 密碼子改變和缺失; 14. 密碼子改變和插入; 15. 同義編碼區(qū); 16. 非同義編碼區(qū); 17. 同義終止區(qū); 18. 終止區(qū); 19. 終止缺失; 20. 其他。

1. Intergenic; 2. Intragenic; 3. Intron; 4. Upstream; 5. Downstream; 6. Splice site acceptor; 7. Splice site donor; 8. Splice site region; 9. Start lost; 10. Frame shift; 11. Condon deletion; 12. Condon insertion; 13. Condon change and condon deletion; 14. Condon change and condon insertion; 15. Synonymous coding; 16. Non synonymous; 17. Synonymous stop; 18. Stop gained; 19. Stop lost; 20. Other.et al.，2013）軟件對可變剪接類型和表達量進行分析，結果表明在9個樣品中5′端和3′端外顯子可變剪切最多（圖2）。

2.2.4 基因結構優(yōu)化 由于處理數(shù)據(jù)的軟件不同或數(shù)據(jù)本身的局限性，本文對中國種茶樹基因組的261個基因結構進行了優(yōu)化并列出了部分優(yōu)化的基因（表4）。

2.3 新基因

2.3.1 新基因功能注釋 我們對利用StringTie拼接出來的轉錄本與中國種茶樹基因組的注釋信息進行比較，發(fā)掘出了13 531個新基因，并將它們與Swiss-Prot，COG，Pfam，KEGG，GO和NR數(shù)據(jù)庫進行BLAST，獲得新基因的注釋信息（表5）。

2.3.2 基因表達量 本文通過最大流量算法，采用FPKM （Fragments Per Kilobase of transcript per1. 可變5′或3′端剪切; 2. 單內含子滯留; 3. 多內含子滯留; 4. 多外顯子跳躍; 5. 單外顯子跳躍; 6. 第一個外顯子可變剪切; 7. 最后一個外顯子可變剪切; 8. 可變5′或3′端剪切（模糊邊界）; 9. 單內含子滯留（模糊邊界）; 10. 多內含子滯留（模糊邊界）; 11. 多外顯子跳躍（模糊邊界）; 12. 單外顯子跳躍（模糊邊界）。

1. Ae; 2. Ir; 3. Mir; 4. Mskip; 5. Skip; 6. Tss; 7. Tts; 8. Xae; 9. Xir; 10. Xmir; 11. Xmskip; 12. Xskip.

福鼎大白種茶根莖葉三個部位的基因表達量的密度分布對比（圖3）表明，能夠測序到的編碼蛋白的基因FPKM值跨越10-2到104六個數(shù)量級。

2.4 差異表達分析

2.4.1 差異表達篩選 利用Deseq（Wang et al.，2010）進行樣品間的差異表達分析，并且將Fold Change≥2，F(xiàn)DR<0.01作為篩選標準，該文列出了部分的差異表達基因并且統(tǒng)計了根莖葉兩兩樣品間的差異基因的數(shù)目（表6，表7）。從葉和莖、葉和根、莖和根的差異表達火山圖（圖4）中可以看出，差異基因的統(tǒng)計學顯著性以及上調基因和下調基因，紅點代表上調，綠點代表下調。此外，我們還計算了三個組織差異表達基因的聚類熱圖（圖5）。

2.4.2 差異基因功能注釋和富集分析 對差異基因進行COG，GO，KEGG，KOG，NR，Pfam，Swiss-Prot和eggNOG功能注釋（表8），并且繪制了葉和莖、 葉和根、 莖和根差異基因的GO功能分類圖（圖6）。應用超幾何檢驗，對差異表達基因進行KEGG通路富集分析，統(tǒng)計了顯著性q值最小的前20個通路（圖7），圖中每個圓圈代表一個pathway通路，圓圈大小代表通路總富集的基因數(shù)目，圓圈越大，表示基因越多。圓圈顏色代表q值，q值越小，代表差異基因在該通路中的富集性越可靠，所以說該圖越靠近右下角的圓圈，參考價值越大。

2.5 DEU分析

本文應用DEXSeq（Anders et al.，2012）進行外顯子水平的差異分析，即DEU（differential exon usage）分析，設置FDR<0.01，該文列出了葉和莖的部分DEU（表9）。

3 討論

福鼎大白種茶樹是異花授粉的經濟作物，基因組雜合度高，具有優(yōu)良的發(fā)芽率、抗寒性和耐旱性（劉本英等，2008）。本研究對都勻地區(qū)福鼎大白種茶根莖葉進行了轉錄組測序，Q30堿基百分比不小于93.22%，樣品與參考基因組的比對效率在87.83%～91.14%之間。李明璽等（2018）對靜安白茶芽和葉的轉錄組研究中，Q30值都為87.32%，Liu et al.（2017）對瑞雪的轉錄組研究中，Q20為97%左右，Wei et al.（2018）對龍井43的轉錄組研究中，Q30為90%左右，因此，福鼎大白種1. 代謝過程; 2. 細胞過程; 3. 單組織過程; 4. 生物調節(jié); 5. 應激反應; 6. 定位; 7. 細胞成分組織或生物合成; 8. 發(fā)育過程; 9. 多細胞生物過程; 10. 信號; 11. 生殖; 12. 生殖過程; 13. 多組織過程; 14. 生長; 15. 解毒作用; 16. 免疫系統(tǒng)過程; 17. 生物附著; 18. 生物相; 19. 節(jié)律性過程; 20. 細胞殺傷性; 21. 運動; 22. 細胞; 23. 細胞組件; 24. 細胞器; 25. 細胞膜; 26. 細胞膜組件; 27. 細胞器組件; 28. 大分子復合物; 29. 胞外區(qū); 30. 細胞膜內控; 31. 細胞連接; 32. 共質體; 33. 超分子復合物; 34. 擬核; 35. 病毒粒子; 36. 病毒粒子組件; 37. 胞外區(qū)組件; 38. 催化活性; 39. 結合; 40. 轉運活性; 41. 結構分子活性; 42. 核酸結合轉錄因子活性; 43. 電子載體活性; 44. 信號轉導活性; 45. 分子功能調節(jié)因子; 46. 抗氧化活性; 47. 分子轉導活性; 48. 營養(yǎng)庫活性; 49. 轉錄因子活性，蛋白結合; 50. 蛋白標簽; 51. 金屬伴侶活性; 52. 翻譯調控活性。

1. Metabolic process; 2. Cellular process; 3. Single-organism process; 4. Biological regulation; 5. Response to stimulus; 6. Localization; 7. Cellular component organization or biogenesis; 8. Developmental process; 9. Multicellular organismal process; 10. Signaling; 11. Reproduction; 12. Reproductive process; 13. Multi-organism process; 14. Growth; 15. Detoxification; 16. Immune system process; 17. Biological adhesion; 18. Biological phase; 19. Rhythmic process; 20. Cell killing; 21. Locomotion; 22. Cell; 23. Cell part; 24. Organelle; 25. Membrane; 26. Membrane part; 27. Organelle part; 28. Macromolecular; 29. Extracellular region; 30. Membrane-enclosed lumen; 31. Cell junction; 32. Symplast; 33. Supramolecular complex; 34. Nucleoid; 35. Virion; 36. Virion part; 37. Extracellular region part; 38. Catalytic activity; 39. Binding; 40. Transporter activity; 41. Structural molecule activity; 42. Nucleic acid binding transcription factor activity; 43. Electron carrier activity; 44. Signal transducer activity; 45. Molecular function regulator; 46. Antioxidant activity; 47. Molecular transducer activity; 48. Nutrient reservoir activity; 49. Transcription factor activity，protein binding; 50. Protein tag; 51. Metallochaperone activity; 52. Translation regulator activity.

1. 有機含硒化合物代謝; 2. 類胡蘿卜素生物合成; 3. 類黃酮生物合成; 4. 糖尿病并發(fā)癥中的年齡信號通路; 5. 淀粉和蔗糖代謝; 6. 碳代謝; 7. 葉酸碳池; 8. 半乳糖代謝; 9. 角質、亞氨酸和蠟的生物合成; 10. 單萜類生物合成; 11. 甘氨酸、絲氨酸和蘇氨酸代謝; 12. 光合生物固碳作用; 13. 乙醛酸和二羧酸代謝; 14. 氮代謝; 15. 氰基氨基酸代謝; 16. 植物激素信號轉導; 17. 卟啉和葉綠素代謝; 18. 苯丙素的生物合成; 19. 光合作用-天線蛋白; 20. 光合作用; 21. 抗壞血酸和醛酸代謝; 22. 硫代謝; 23. 雙醌和其它萜烯醌生物合成; 24. 谷胱甘肽代謝; 25. 萜類骨架生物合成; 26. 磷酸戊糖途徑; 27. 脂肪酸降解; 28. 糖酵解/糖新生; 29. 亞麻酸代謝; 30. 二萜生物合成; 31. 倍半萜和三萜生物合成。

1. Selenocompound metabolism; 2. Carotenoid biosynthesis; 3. Flavonoid biosynthesis; 4. AGE-RAGE signalling pathway in diabetic complications; 5. Starch and sucrose metabolism; 6. Carbon metabolism; 7. One carbon pool by folate; 8. Galactose metabolism; 9. Cutin，suberine and wax biosynthesis; 10. Monoterpenoid biosynthesis; 11. Glycine，serine and threonine metabolism; 12. Carbon fixation in photosynthetic organism; 13. Glyoxylate and dicarboxylate metabolism; 14. Nitrogen metabolism; 15. Cyanoamino acid metabolism; 16. Plant hormone signal transduction; 17. Porphyrin and chlorophyll metabolism; 18. Phenylpropanoid biosynthesis; 19. Photosynthesis-antenna proteins; 20. Photosynthesis; 21. Ascorbate and aldarate metabolism; 22. Sulfur metabolism; 23. Biquinone and other terpenoid-quinone biosynthesis; 24. Glutathione metabolism; 25. Terpenoid backbone biosynthesis; 26. Pentose phosphate pathway; 27. Fatty acid degradation; 28. Glycolysis/Gluconeogenesis; 29. Alpha-Linolenic acid metabolism; 30. Diterpenoid biosynthesis; 31. Sesquiterpenoid and triterpenoid biosynthesis.

茶文庫構建成功且測序質量良好，可進行下一步分析。SNP位點數(shù)目為489 642～638 706，轉換類型比例為65.46%～66.33%，顛換類型為33.67%～34.54%，雜合型SNP位點比例為34.70%～51.38%。對261個基因進行了結構優(yōu)化，發(fā)現(xiàn)了13 531個新基因，其中10 244個新基因被注釋了，注釋到COG數(shù)據(jù)庫中有2 363個，注釋到GO的有6 218個，注釋到KEGG的有3 948個，注釋到KOG的有5 833個，注釋到Pfam的有5 383個，注釋到Swiss-Prot的有6 900個，注釋到eggNOG的有9 243個，注釋到NR的有10 107個。

福鼎大白種茶根莖葉中有12 595個基因參與代謝過程，12 507個基因參與催化活性，8 176個基因參與細胞組成成分，是基因數(shù)最多的三個生物過程。陳琳波等（2015）對“紫鵑”茶樹的轉錄組分析中，Unigene的GO分析表明，代謝過程、細胞、細胞部分、催化活性等富集程度較高，與本研究結果相一致。三個組織中有 1 615個基因參與發(fā)育過程（278，405，191），即葉和莖參與此過程的有278個差異基因，葉和根中有405個差異基因，莖和根中有191個。根莖葉有1 401個基因參與運輸活動（197，356，231），288個基因參與生長過程（47，65，20），262個基因參與信號轉導活動（44，70，30），200個基因參與免疫反應過程（35，56，36）。

都勻地區(qū)茶獨特的風味是由葉片中的茶氨酸、萜烯類和類黃酮的代謝產物賦予的。Wu et al.（2013）從茶樹葉片的轉錄組研究中得出大多數(shù)基因參與編碼合成黃酮類、咖啡堿和茶氨酸等次生代謝途徑的重要酶。該文通過KEGG功能顯著性富集分析，發(fā)現(xiàn)葉和莖、葉和根、莖和根差異基因參與的部分極其顯著的通路有類胡蘿卜素生物合成，泛醌和其他萜烯類醌的生物合成，單萜生物合成，甘氨酸、絲氨酸和蘇氨酸代謝，乙醛酸和二羧酸代謝，植物激素信號轉導，卟啉和葉綠素代謝，苯丙素的生物合成，光合作用天線蛋白和光合作用，類黃酮生物合成，萜類化合物的生物合成，亞麻酸代謝等顯著富集。葉和莖的差異表達基因在KEGG通路中，有90個基因參與碳代謝（8.65%，21，69），占所有差異基因的8.65%，其中有21個差異基因在莖中表達上調，69個差異基因在葉中表達上調。79個差異基因參與植物激素信號轉導（7.6%，55，24），75個差異基因參與苯丙素的生物合成 （7.21%，49，26），70個差異基因參與淀粉和蔗糖代謝（6.73%，46，24，），64個差異基因參與氨基酸生物合成（6.15%，22，42），43個差異基因參與植物病原體互作通路（4.13%，30，13），19個差異基因參與類黃酮的生物合成（1.83%，13，6），20個差異基因參與萜類化合物生物合成（1.92%，4，16），24個差異基因參與卟啉和葉綠素代謝（2.31%，1，23），44個差異基因參與光合作用（4.23%，0，44）。以上研究表明都勻地區(qū)福鼎大白種茶獨特的香味和各種健康功能主要歸因于萜類物質合成、光合作用、 類黃酮的生物合成等生物過程，類黃酮主要賦予茶的苦味，二萜類化合物的水解產物具有令人愉悅的香味和特色風味。為進一步研究都勻地區(qū)福鼎大白種茶生長發(fā)育和遺傳育種提供了分子基礎，并為揭示與類黃酮相關代謝通路有關基因提供了重要理論依據(jù)。