海菜花轉錄組測序及生物信息學分析

李健玲，秦 波，黃 欣，蔣日紅，孫 苗，梁圣華，黃耀恒，韋廣綏

（1.廣西壯族自治區林業科學研究院，廣西特色經濟林培育與利用重點實驗室，廣西 南寧530002；2.北京林業大學，國家花卉工程技術研究中心，北京 100083；3.廣西國有高峰林場，廣西 南寧 530025）

海菜花（Ottelia acuminata）屬水鱉科水車前屬多年生沉水植物，為中國特有種，國家二級保護植物。該植物對水質要求較高，是一種環境指示性植物，主要分布區域為廣東、海南、廣西、四川、貴州、云南等地[1]。海菜花是一種典型的高鉀低鈉型蔬菜，含有鈣、鐵、蛋白質、抗壞血酸、酚類等豐富的營養物質，具有較高的食用價值；其中酚類物質是天然的抗氧化劑，對DNA 損傷有良好的修復作用。另一方面，海菜花還能富集重金屬元素鉛，具有一定的生態修復價值[2-3]。目前，關于海菜花的研究多集中在栽培[4-5]、進化和遺傳學分析[6-7]、化學成分分析[3]等方面，其分子水平的研究僅有葉綠體基因組的報道[8]。

分子生物學技術在植物領域尤其是農作物方面應用的較早且深入，但是在水生植物及藥用植物上的研究報道仍然有限。高通量轉錄組作為初步掌握代謝通路和生物合成基因信息的一個手段，在許多物種中已經得到了廣泛的應用。為了彌補海菜花轉錄組上的空白，研究對海菜花葉片進行轉錄組測序，挖掘其重要基因和調控信息，以期為海菜花功能基因和遺傳多樣性分析及分子育種奠定基礎。

1 材料與方法

1.1 試驗材料

采集新鮮的海菜花葉片，并立即投入液氮中保存。

1.2 海菜花cDNA 文庫構建及測序

使用Trizol 法對海菜花的RNA 進行提取。利用Illumina TruSeqTM RNA sample prep Kit（Illumina，美國）方法構建海菜花的RNA 文庫，首先使用帶有Oligo（dT）的磁珠對帶有polyA 尾巴的mRNA 進行富集，并用超聲波將mRNA 進行片段化打斷。加入隨機寡核苷酸為引物，在M-MuLV 逆轉錄酶體系下反轉合成第一條鏈的cDNA，之后用RNaseH 清除反應中的RNA，在DNA polymerase I 體系下，以dNTPs 為原料進行cDNA 第二條鏈的合成。對合成后的雙鏈cDNA 進行純化、末端修復、加A 尾以及連接測序接頭，利用AMPure XP beads 對處理過的雙鏈cDNA 進行篩選，挑出200 bp 左右的片段進行擴增，之后對得到的PCR 產物進行純化，最終產物即為文庫。

利用瓊脂糖凝膠電泳、NanoPhotometer spectrophotometer、Qubit2.0 Fluorometer、Agilent 2100 bioanalyzer 對RNA 的質量進行檢驗，樣品合格后進行測序。在測序的flow cell 中加入不同標記的4 種dNTP 以及DNA 聚合酶和接頭引物，當測序鏈在延伸互補的時候，每個dNTP 會發出相應的熒光，從而被測序儀捕捉獲得序列信息。

海菜花的RNA 提取、文庫構建及測序均由上海凌恩生物科技有限公司完成。

1.3 轉錄組組裝及分析

使用Trimmomatic 對原始測序數據進行質控，Trinity 進行轉錄組拼接，并使用Salmon 進行拼接數據的比對，之后使用Diamond 對測序數據進行NR數據庫、KEGG 數據庫、Swiss-Port 數據庫的注釋，利用eggNOG-mapper 進行eggNOG 數據庫和GO 數據庫的注釋。

2 結果與分析

2.1 海菜花轉錄組測序及組裝

對海菜花的葉片進行轉錄組測序，共檢測到26 658 026 條reads，共計3 998 703 900 個堿基。對原始數據進行質控進一步去除reads 中的接頭序列、質量較低的堿基以及短序列后，得到clean data 的reads 有 26 497 174 條，共計3 961 376 102 個堿基，其中Q20%高達98.90%，Q30%高達96.19%，GC含量占比為50.57%，表明海菜花轉錄組數據質量較高，可以開展后續的數據分析。

由于海菜花沒有參考基因組，因此利用Trinity（2.8.6）對clean data 進行從頭組裝，得到組裝后的unigene 有27 040 條，長度為27 761 688 bp，最長的unigene 有12 094 bp，N50 為1 340 bp，平均unigene 為1 026.69 bp，其中N50 高于平均長度說明組裝效果較好。同時，將測序獲得的高質量序列與拼接序列進行比對，比對率達到86.89%，也表明拼接效果較好。在得到的unigene 中，401～600 bp 的序列長度占比最大，達到7 134 條（26.38%），序列長度分布見圖1。

圖1 海菜花unigene 長度分布

2.2 海菜花轉錄組注釋

將獲得的27 040 條unigene 分別與不同數據庫進行比對，結果如圖2 所示。有18 910 條unigene分別被NR、GO、COG、KEGG、SWISS 這5 大數據庫注釋到，占unigene 總數的69.93%。其中，NR數據庫比對到的unigene 最多，有18 874 條，占總數的69.80%；其次是SWISS 數據庫，比對到14 020條unigene，占總數的51.85%；第三是COG 數據庫，比對到13 531 條unigene，占總數的50.04%；KEGG 數據庫比對到的unigene 數量最少，為7 552條，占總數的27.93%；在所有數據庫中均比對到的unigene 有4 819 條，占總數的17.82%，推測為新基因。

圖2 海菜花unigene 在各數據庫的注釋結果

2.3 NR 數據庫注釋

NR數據庫的功能注釋能夠了解海菜花轉錄組中序列與哪些物種有同源分布。由圖3 可知，海菜花與搖蚊（Clunio marinus）比對到的同源序列最多，為2 622 條，占所有比對到序列的13.89%；其次是油棕（Elaeis guineensis），達到2 062 條，占10.93%；海棗（Phoenix dactylifera）中比對到1 646條，占8.72%；另外，芋（Colocasia esculenta）、大葉藻（Zostera marina）、蓮（Nelumbo nucifera）、鳳梨（Ananas comosus）、小果野蕉（Musa acuminatasubsp.malaccensis）、沉水樟變型（Cinnamomum micranthumf.kanehirae）、石刁柏（Asparagus officinalis）中 分 別 比 對 到1 435、947、750、637、570、570、544 條，分別占比7.60%、5.02%、3.97%、3.38%、3.02%、3.02%、2.88%；而其他物種中比對到的序列數均低于500 條，總數達7 091 條，占比為37.57%。

圖3 海菜花unigene 在NR 數據庫中比對到的同源序列排名前10 的物種

2.4 COG 數據庫注釋

由圖4 可知，海菜花轉錄組數據在COG 數據庫中共比對到unigene 13 531 條，但功能未知的基因占多數，排在KOG 預測的第1 位，表明海菜花中有許多序列還未得到有效的驗證和研究；其中，轉錄后修飾、蛋白質轉換、伴侶蛋白（O：Posttranscriptional modification，protein turnover，chaperones），信號傳導機制（T：Signal transduction mechanisms），翻譯、核糖體結構與生物發生（J：Translation，ribosome structure and biogenesis），RNA轉錄和修飾（A：RNA processing and modification）分別排在KOG 預測的第2～5 位，代表了海菜花生長發育和信號轉導中基本過程活動比較活躍。

圖4 海菜花葉片轉錄組COG 功能注釋

2.5 GO 注釋

在GO 注釋中，海菜花共注釋到10 505 條unigene，分為生物過程（biological process）、細胞組分（cellular component）和分子功能（molecular function）3 個大類和43 個功能亞類（圖5）。生物過程注釋到的unigene 最多，亞類中細胞進程和代謝過程中unigene 數量最多，而碳利用和生物礦化所占比例最少。在分子功能中，催化活性和結合比對到的unigene 數量最多，蛋白標簽、運貨受體活性、養分庫活性和小分子傳感器活性比對到的數量較少。在細胞組分中包括了4 個亞類，最高的為細胞解刨學實體，最少的為其他器官部分。

圖5 海菜花葉片轉錄組GO 注釋結果

2.6 KEGG 通路分析

為了研究基因具體的功能，對海菜花轉錄組中的unigene 進行KEGG 通路分析，共分為5 大類和20 個亞類（圖6）。第1 大類代謝過程中比對到的unigene 數量最多，包括11 個亞類，全局和概述圖譜比對到的數量最多，達4 678 條unigene，其次是碳水化合物代謝，比對上1 047 條unigene，氨基酸代謝、能量代謝、酯類代謝分別比對到587、498、486 條，表明在細胞過程中關于代謝的基因較多。第2 大類是遺傳信息處理，包括5 個亞類，翻譯過程中比對到的unigene 最多，達到784 條，其次是折疊、分類和降解，比對到662 條，病毒信息中比對到的數量最少，為46 條。第3 大類為環境信息處理，包括2 個亞類，分別為信號轉導（比對到309條）和膜轉運（比對到25 條）。第4 大類為細胞過程，包括2 個亞類，為運輸和分解代謝（比對到427 條）和細胞運動（比對到84 條）。最后1 類為有機系統，僅包含1 個亞類，為環境適應，比對到200 條unigene。

圖6 海菜花葉片轉錄組KEGG 注釋分類統計

2.7 SSR 分析

在海菜花轉錄組中共檢索到4 217 個SSR 位點，包含了6 種核苷酸重復類型，具體情況如圖7 所示。單核苷酸重復數量最多，為1 518 個，占比為36.00%，A/T 重復比例最高；其次為三核苷酸重復，有1 437 個，占比為34.08%，且類型較多；雙核苷酸重復為1 214 個，占比為28.79%，以AG/CT 類型比例最高；四核苷酸位點重復有34 個，占比為0.81%；六核苷酸位點重復有9 個，占比為0.21%；五核苷酸位點重復數最少，僅5 個，占比0.11%。

圖7 海菜花葉片轉錄組SSR 分布情況

3 討 論

海菜花在景觀營造、食用和藥用方面有很高的開發利用價值。目前，對水生植物轉錄組的研究種類較少，包括圓海鏈藻（Thalassiosira rotula）[9]、蓮（Nelumbo nucifera）[10-11]、齒葉睡蓮（Nymphaea lotus）[12]、東方澤瀉（Alisma orientale）[13]、芡實（Euryale feroxSalisb.ex Konig et Sims）[14]等，通過轉錄組測序技術對其次生代謝物的合成、花開放規律、脅迫響應和遺傳標記等進行了深入研究。海菜花在生態適應性上對環境的要求較高，且有不同的生態適應類群，次生代謝物豐富，對其進行分子層面的深入研究具有重要意義。

該研究在高通量轉錄組測序技術的支持下，對海菜花葉片的轉錄組信息進行解讀，共檢測到unigene 27 040 條，其中有18 910 條unigene 在數據庫中比對上了，獲得了海菜花中大量的基因信息。在NR 數據庫得到了最多的基因注釋，共18 874 個，與植物同源基因比對最多的是油棕，另外包括芋、大葉藻、蓮和鳳梨等，這些植物類型均屬于濕生環境，表明海菜花可能在轉錄調控和生態適應方面與濕生植物有相似的基因調控類型。在KOG 數據庫中，海菜花與已知基因比對上最多的基因涉及轉錄后修飾、蛋白質轉換、伴侶蛋白等方面。對基因通路進行KEGG 分析，共比對到7 552 個unigene，以細胞過程比對到的unigene 最多，且與代謝途徑密切相關；另外，在GO注釋中，生物過程注釋到的unigene最多，表明在海菜花葉片中新陳代謝活動較為旺盛。SSR分析中單核苷酸、雙核苷酸和三核苷酸重復是主要類型，占比達98.87%，這些位點的發現對海菜花的分子標記開發和利用具有重要價值，能為植物鑒定、遺傳育種和藥理研究提供直接的證據。

通過對海菜花葉片進行高通量轉錄組測序，獲得了大量的基因序列信息和注釋信息，為海菜花的功能基因和遺傳多樣性分析及分子育種奠定了基礎，也為海菜花的生態適應性、食用和藥用價值的進一步開發利用提供了依據。