厚葉木蓮(Manglietia pachyphylla)基因組草圖*

甘新軍,賓 粵,陳煥錦,朱韋光,熊露橋,余恩萍,4,王崢峰**,徐鳳霞,曹洪麟

(1.廣東從化陳禾洞省級自然保護區管理處,廣東廣州 510950;2.中國科學院華南植物園,廣東省應用植物學重點實驗室,中國科學院退化生態系統植被恢復與管理重點實驗室,廣東廣州 510650;3.華南國家植物園,廣東廣州 510650;4.中國科學院大學,北京 100049;5.中國科學院華南植物園,中國科學院植物資源保護與可持續利用重點實驗室,廣東廣州 510650)

木蘭科(Magnoliaceae)是最原始的被子植物,其屬種豐富,類型多樣,是研究被子植物起源和進化的重要類群[1-4]。我國是木蘭科起源地和避難所,為木蘭科植物分布中心,擁有很多古老、孑遺和特有物種[2,5-7]。我國的木蘭科物種現主要分布于我國熱帶亞熱帶地區,如云南省、廣西壯族自治區、貴州省和湖南省。木蘭科中的很多物種樹形優美,葉型秀麗、色澤鮮艷,花形態各異、花色豐富明艷且高貴典雅,觀賞性強,是優良的綠化樹種[8-11]。木蘭科物種的枝、葉、花含有揮發性有機物,作為綠化樹種可以凈化空氣,促進人們身心健康[12-15],而且這些揮發性有機物與木蘭科產生的其他次生代謝物如生物堿,常被用作中藥[16-18]。另外,木蘭科的一些種類具有較強的光合能力、固碳能力以及土壤改良作用,可用于人工林改造,實現林業提質增效[19,20];而且一些木蘭科植物樹干通直,木材材質均勻、細密,是很好的用材樹種[ 8,21,22]。同時,木蘭科植物最早記載于我國秦漢時代,歷經千年,寄托了人們對美好生活的追求和向往,富有深厚的文化沉淀[23]。

厚葉木蓮(Manglietiapachyphylla)為木蘭科木蓮屬(Manglietia)的木本植物,為中國特有種,目前零星分布于我國廣東省和廣西壯族自治區海拔500 m以上的常綠闊葉林中[24,25]。厚葉木蓮的顯著特征是有光澤的革質厚葉,可作為園林綠化和用材物種[24]。厚葉木蓮種群少,個體數量小,結實率低,種子易被動物啃食,自然更新差[25,26],現被列為國家二級重點保護野生植物[24]。目前國內外已開展厚葉木蓮群落學[25,26]、花粉形態[27]和光響應生理[28]等方面的研究,但還未有關于厚葉木蓮的遺傳多樣性及其基因組的研究報道?；蚪M及基于基因組開展的物種遺傳多樣性研究,可以揭示物種進化過程,從而了解物種的適應性并進一步應用于品種改良[29]。為此,本研究采用二代和三代高通量測序手段,對厚葉木蓮基因組進行測序,組裝其基因組草圖,為今后更好地開展其進化、遺傳多樣性研究提供參考。

1 材料與方法

1.1 材料

在廣東省廣州市從化區陳禾洞省級自然保護區內選擇1株厚葉木蓮成年大樹(生長點地理位置為113°55′31.84″E,23°45′1.25″N),其胸高直徑(Diameter at Breast Height,DBH)為22.4 cm,采集其無蟲咬和病斑的樹葉3片。采集的樹葉用剪刀剪碎后,用錫箔紙包裹后立刻投入液氮罐中,之后將樣品送至武漢未來組生物技術有限公司進行高通量測序。

1.2 方法

測序包括3個方面的內容:一是采用Nanopore PromethION測序平臺對厚葉木蓮進行三代基因組測序,二是采用MGI DNBSEQ-T7 測序平臺對厚葉木蓮進行二代基因組測序,三是采用MGI DNBSEQ-T7 測序平臺對厚葉木蓮進行轉錄組測序,具體實驗流程參考Wang等[30]的研究。針對所得的測序數據,利用不同程序開展數據處理和基因組組裝、分析。在數據處理過程中主要使用程序的默認參數,如在程序運行過程中對默認參數進行改動,則會在文中具體說明。

1.2.1 測序數據前處理

本研究利用Sickle v1.33 (https://github.com/najoshi/sickle)對厚葉木蓮二代基因組測序數據進行過濾,去除測序數據中堿基質量小于30、片段長度小于80 bp的測序數據。過濾后的二代測序數據用 RECKONER v1.1[31]進行糾錯。對于三代基因組測序數據,由于低質量數據在交付前測序公司已過濾,無需再過濾,本研究只利用Porchop 0.2.4 (https://github.com/rrwick/Porechop)對厚葉木蓮三代基因組測序數據進行接頭過濾。

1.2.2 基因組大小預測與組裝

針對厚葉木蓮二代測序數據,利用GenomeScope 2.0[32]進行基因組大小預測(參數“-k 21”)。針對厚葉木蓮三代測序數據,本研究選擇大于10 kb測序序列,利用NextDenovo 2.3.1(https://github.com/Nextomics/NextDenovo)進行基因組組裝,組裝的基因組利用Pseudohaploid (https://github.com/schatzlab/pseudohaploid)和Purge_Dups v1.2.6[33]去除冗余序列(如雜合導致的拼接序列),之后利用racon v1.5.0[34],hapo-G v1.3.2[35]和polypolish v0.5.0[36]進行組裝序列糾錯。針對組裝完成的基因組,利用BUSCO (Benchmarking Universal Single-Copy Orthologs)v5.4.6[37]對照“eudicots_odb10.2020-09-10”和“embryophyta_odb10.2020-09-10”兩個單拷貝基因庫進行組裝完整性的評估。

1.2.3 重復序列和基因的預測、注釋

本研究利用EDTA (Extensive De-novo TE Annotator)v2.1.0[38]和RED (REpeat Detector)v2.0[39]預測厚葉木蓮基因組中的重復序列。利用BEDTools v2.29.2[40]中的“merge”命令將兩個重復序列預測結果合并,再利用“maskfasta”命令將預測的重復序列屏蔽掉。針對屏蔽了重復序列的基因組,首先利用BRAKER2[41],同時結合轉錄組數據和9個物種的蛋白質序列(表1)預測厚葉木蓮基因組組裝序列中的基因;然后將結果輸入funannotate pipeline v1.8.13(https://github.com/nextgenusfs/funannotate)中,同樣結合轉錄組數據和9個物種的蛋白質序列,進一步預測厚葉木蓮基因組組裝序列中的基因,這一過程包括3個步驟和命令:“funannotate train”“funannotate predict”和“funannotate update”,在“predict”步驟中使用的參數為“-max_intronlen 100,000 -busco_db embryophyta -organism other”。

表1 預測厚葉木蓮基因組基因的參考物種

基因預測結束后,利用8個不同的蛋白質注釋平臺對預測的基因進行功能分析,包括:dbCAN (DataBase for automated Carbohydrate-active enzyme ANnotation)v10.0[42],eggNOG-mapper (Evolutionary genealogy of genes:Non-supervised Orthologous Groups-mapper)v5.0.2[43],GO (Gene Ontology)[44,45],KEGG (Kyoto Encyclopedia of Genes and Genomes)[46],InterPro (The Integrated Resource of Protein Domains and Functional Sites)v5.60-92.0[47],MEROPS v12.0[48],Pfam (The protein families database)v35.0[49]和SignalP 5.0b[50]。

考慮到注釋的基因中有可變剪切的狀況,在進行比較基因組的研究中去除了各物種基因中的可變剪切產生的基因,只保留其中最長的基因序列進行分析。

1.2.4 基因家族和系統發育分析

利用OrthoFinder 3.0.0[51,52]并結合其他9個物種(表1)的蛋白質序列進行基因家族分析。在OrthoFinder分析過程中,程序自動分析獲得物種間的單拷貝基因,并用這些單拷貝基因進行系統發育樹構建。利用構建的系統發育樹,使用MCMCTree[53]進行物種間分化時間的估算。在這一分析過程中需要參考已有物種間的分化時間,本研究從http://timetree.org/獲得這些信息,并在表2中列出。在得到有分化時間的系統發育樹后,再利用cafe (Computational Analysis of gene Family Evolution)v5[54]進行基因家族的擴張和收縮分析;對其中顯著擴張和收縮的基因家族,則利用 TBtools v1.115[55]進行GO和KEGG的富集分析。

表2 從http://timetree.org/獲得的物種對間分化時間

1.2.5 基因重復(Gene duplications)

利用wgd v1.1.2[56]進行基因組的全基因組重復(Whole genome duplication)事件分析,該過程是在基因組中進行基因間的同源性分析,找到兩兩同源基因,并對兩兩同源基因進行同義突變率(synonymous substitution rate,Ks)計算,之后查看Ks值的密度分布圖,其峰值出現的地方提示有全基因重復事件發生;進一步利用 Doubletrouble v0.99.1 (https://github.com/almeidasilvaf/doubletrouble)開展全基因組重復基因、串聯重復(Tandem duplications)基因、近端重復(Proximal duplications)基因、轉座重復(Transposed duplications)基因、散在重復(Dispersed duplications)基因的分析[55]。其中,串聯重復基因是指彼此連續或中間間隔不超過5個其他基因的近緣基因;近端重復基因是連續分布、中間間隔不超過10個其他基因的近緣基因;轉置重復基因是由轉座子介導的重復基因;散在重復基因是隨機分布、彼此間不靠近的重復基因[57]。利用TBtools分別對全基因組重復基因、串聯重復基因、近端重復基因進行GO和KEGG富集分析。

2 結果與分析

2.1 厚葉木蓮基因組測序與組裝

本研究中三代基因組測序共獲得大約118.1 Gb測序數據,二代基因組測序獲得約264.1 Gb測序數據,轉錄組測序獲得約31.4 Gb測序數據?；蚪M和轉錄組原始測序數據已上傳至GenBank,三代基因組測序數據序列號為SRR24423593、SRR24423594、SRR24423595;二代基因組測序數據序列號為SRR24390471、SRR24390472、SRR24390473;轉錄組測序數據序列號為SRR24415003。利用GenomeScope 2.0分析厚葉木蓮基因組大小為1 969 269 649 bp。

本研究利用NextDenovo 2.3.1進行厚葉木蓮基因組組裝,得到基因組組裝大小為2 350 821 062 bp,包含1 436個拼接序列(contig),N50(將組裝的序列按照長度由大到小進行累加,當累加到某個序列時,累加的值為基因組50%的長度時,此序列的長度即為N50)為6 839 193 bp,最長拼接序列為26 073 671 bp,最短拼接序列為16 675 bp,平均為1 637 062.0 bp。去除冗余序列和糾錯后,最終的基因組組裝大小為2 092 298 891 bp,包含676個拼接序列,N50為7 961 115 bp,最長拼接序列為26 180 362 bp,最短拼接序列為27 281 bp,平均為3 095 117 bp,組裝的其他統計信息見表3。組裝的厚葉木蓮基因組數據上傳至GenBank,序列號為JASAUF000000000。

表3 厚葉木蓮基因組組裝結果

BUSCO對組裝的基因組完整性的評估顯示,比對“eudicots_odb10.2020-09-10”單拷貝基因庫,全部2 326個BUSCO單拷貝基因,96.6% BUSCO單拷貝基因在厚葉木蓮基因組中完整匹配。其中,能完整匹配到BUSCO單拷貝基因庫而且在厚葉木蓮基因組中也是單拷貝的基因占90.8%,能完整匹配到BUSCO單拷貝基因庫但在厚葉木蓮基因組中為多拷貝的基因占5.8%;不完整匹配到BUSCO單拷貝基因庫的基因有27個,占1.2%;有54個BUSCO單拷貝基因未在厚葉木蓮基因組的基因中匹配到,占2.3%。比對“embryophyta_odb10.2020-09-10”單拷貝基因庫,全部1 614個BUSCO單拷貝基因中,98.8% BUSCO單拷貝基因在厚葉木蓮基因組中完整匹配,能完整匹配到BUSCO單拷貝基因庫而且在厚葉木蓮基因組中也是單拷貝的基因占95.3%,能完整匹配到BUSCO單拷貝基因庫但在厚葉木蓮基因組中為多拷貝的基因占3.5%;不完整匹配到BUSCO單拷貝基因庫的基因有12個,占0.7%;有7個BUSCO單拷貝基因未在厚葉木蓮基因組的基因中匹配到,占0.4%。

2.2 厚葉木蓮基因組注釋

在厚葉木蓮基因組中,RED和EDTA分別檢測到1 342 604 397 bp (64.16%)和1 447 436 839 bp(69.17%)重復序列,兩者合并后得到1 601 108 919 bp的重復序列,占基因組的76.52%。EDTA檢測結果表明,厚葉木蓮基因組中重復序列最多的是Gypsy類的長末端重復序列(Long terminal repeat),有475 508 695 bp,占基因組的22.73%(表4)。

表4 EDTA檢測的重復序列

本研究對厚葉木蓮基因組進行基因預測共獲得37 900個基因,這些基因共編碼了41 675種蛋白質。利用數據庫對這些蛋白質序列進行功能注釋,其中32 249個 (77.4%)蛋白質序列注釋到了8個蛋白質注釋平臺數據庫中的一個,具體如下:17 886個蛋白質序列注釋到GO數據庫,21 430個蛋白質序列注釋到InterPro數據庫,31 146個蛋白質序列注釋到eggNOG-mapper數據庫,23 203個蛋白質序列注釋到Pfam數據庫,1 160個蛋白質序列注釋到dbCAN數據庫,1 016個蛋白質序列注釋到MEROPS數據庫,3 078個蛋白質序列注釋到SignalP數據庫,15 752個蛋白質序列注釋到KEGG數據庫?；蜃⑨屛募焉蟼髦林袊茖W院華南植物園數據存儲庫(https://doi.org/10.57841/casdc.0001214)。

2.3 厚葉木蓮基因家族與基因組進化分析

本研究對包括厚葉木蓮在內的10個物種的蛋白質序列進行了基因家族分析,結果共得到24 616個基因家族。從表5統計結果可以看出,厚葉木蓮基因組中有34 319個基因歸為其中一個基因家族,占總基因的90.6%。對于木蘭科其他兩個物種鵝掌楸和望春玉蘭,它們分配到基因家族的基因比例分別為95.7%和93.4%,均高于厚葉木蓮;但厚葉木蓮基因分布于15 894個基因家族中,占全部基因家族的64.6%;而鵝掌楸和望春玉蘭的基因出現在13 858和14 935個基因家族中,占全部基因家族的56.3%和60.7%,低于厚葉木蓮。24 616個基因家族中有710個厚葉木蓮特有的基因家族,這些家族包含了3 373個厚葉木蓮基因。鵝掌楸和望春玉蘭特有的基因家族分別是437和999個,雖然一個低于厚葉木蓮特有基因家族數目,一個高于厚葉木蓮特有基因家族數目,但是它們所包含的基因數(分別為6 845和4 543)均高于厚葉木蓮,而且這些特有基因數占兩個物種各自所有基因的比例也高于厚葉木蓮,特別是鵝掌楸,為19.4%(表5,圖1)。3種木蘭科植物共同出現的基因家族有345個(圖1)。

圖1 厚葉木蓮與其他物種基因家族交集圖

表5 基因家族及其基因信息統計

對厚葉木蓮特有的基因家族進行GO和KEGG富集分析,結果表明這些基因家族所包含的基因的主要功能與細胞內穩態(GO:0055067、GO:0006885、GO:0030433、GO:0036503)、葡聚糖代謝(Glucan metabolic process)(GO:0044042、GO:0016762、GO:0046527)、原核生物抗性(Prokaryotic defense system)、苯丙氨酸代謝(Phenylalanine metabolism)和轉錄(Transcription)等有關(https://doi.org/10.6084/m9.figshare.23690001.v3,表S1和S2)。

系統發育分析表明厚葉木蓮與望春玉蘭聚在一起(圖2),兩個物種是在約10 500 000年前(95%CI:4 989 810-17 055 600年)從共同的祖先開始分化。對基因家族擴張和收縮分析結果表明,厚葉木蓮基因組中有686個基因家族表現為收縮,1 295個基因家族表現為擴張,其中有136個基因家族表現為顯著收縮,417個基因家族表現為顯著擴張。GO和KEGG富集分析表明,厚葉木蓮基因組中顯著擴張的基因家族與木質部、韌皮部發育(GO:0010088,GO:0010087)、肌動蛋白絲(Actin filament,GO:0030837、GO:0030833、GO:0061572、GO:0051017、GO:0030832、GO:0008064、GO:0030041)、細胞運動(Cell motility)、油菜素甾醇生物合成(Brassinosteroid biosynthesis)、類黃酮生物合成(Isoflavonoid biosynthesis)、維生素B6代謝(Vitamin B6 metabolism)、硫代葡萄糖苷(Glucosinolate biosynthesis)和萜類物質生物合成有關(https://doi.org/10.6084/m9.figshare.23690001.v3,表S3和表S4);厚葉木蓮基因組中顯著收縮的基因家族主要與木質素(Lignin)代謝(GO:0009808、GO:0046274)、類黃酮生物合成,以及二苯乙烯、二芳基庚酸類和姜酚生物合成(Stilbenoid,diarylheptanoid and gingerol biosynthesis)相關(https://doi.org/10.6084/m9.figshare.23690001.v3,表S5和表S6)。

Divergence times are shown below the tree;The "-/+" and the numbers beside the tree nodes and species represent the number of contracted and expanded gene families in Manglietia pachyphylla and other species.

2.4 厚葉木蓮全基因組和基因重復分析

全基因組重復事件分析表明,厚葉木蓮和鵝掌楸、望春玉蘭表現出同樣的Ks峰型(圖3),因此它們近期共同經歷了一次全基因組重復事件。對厚葉木蓮基因組的基因重復研究表明,厚葉木蓮基因組中有4 769個基因與其全基因組重復有關(12.6%),3 317個基因為串聯重復(8.8%),3 124個基因為近端重復(8.2%),136個基因為轉置重復(0.4%),18 522個基因為散在重復(48.9%)。GO和KEGG富集分析表明,厚葉木蓮基因組中全基因組重復基因主要功能與DNA內復制調控(Regulation of DNA endoreduplication)、植物晝夜節律(Circadian rhythm-plant)、鳥嘌呤核苷酸結合蛋白(GTP-binding proteins)、檸檬酸循環[Citrate cycle (TCA cycle)],蛋白酶體(Proteasome)、轉錄因子(Transcription factors)(https://doi.org/10.6084/m9.figshare.23690001.v3,表S7和表S8)相關。GO和KEGG富集分析表明,厚葉木蓮基因組中串聯重復基因與NADH氧化(NADH oxidation,GO:0006116)、寡肽轉運(Oligopeptide transport,GO:0006857)、RNA脫帽(RNA decapping,GO:0110154)、谷胱甘肽代謝過程(Glutathione metabolic process,GO:0006749)、過氧化氫分解代謝過程(Hydrogen peroxide catabolic process,GO:0042743、GO:0042744)、類黃酮生物合成、?；撬岷偷团；撬岽x(Taurine and hypotaurine metabolism)、苯并惡唑嗪酮類化合物生物合成(Benzoxazinoid biosynthesis)、玉米素生物合成(Zeatin biosynthesis)、萜類生物合成(Terpenoid biosynthesis)、苯丙烷類生物合成(Phenylpropanoid biosynthesis)、生物堿合成等相關(https://doi.org/10.6084/m9.figshare.23690001.v3,表S9和表S10)。厚葉木蓮基因組中近端重復基因主要和防御(GO:0043207、GO:0098542、GO:0051707、GO:0009607、GO:0050832、GO:0009605、GO:0006952)、次生代謝(GO:0044550、GO:0019748)、各類生物堿生物合成(Biosynthesis of various alkaloids)、花青素生物合成(Anthocyanin biosynthesis)、聚酮生物合成(Polyketide biosynthesis)、苯并惡唑嗪酮類化合物生物合成、萜類生物合成、硫代葡萄糖苷生物合成和油菜素甾醇生物合成等相關(https://doi.org/10.6084/m9.figshare.23690001.v3,表S11和表S12)。

圖3 全基因組重復分析中同義突變率(Ks)的密度分布

3 討論

3.1 厚葉木蓮的基因組組裝

高通量測序已成為研究物種基因組和遺傳多樣性的主要手段。目前高通量測序分為短片段測序和長片段測序量兩種方式,前者測序平臺主要包括Illumia和MGI測序公司一系列儀器設備,后者測序平臺主要包括Nanopore和PacBio測序公司一系列儀器設備。短片段高通量測序的優勢是測序的準確性較高、數據量大和價格便宜;長片段測序的優勢是測得的序列長度長、完整性好,有利于后續基因組組裝的連續性,如Nanopore測序平臺可測100 kb或更長的序列片段。但采用Nanopore和PacBio測序平臺測得的序列錯誤率較高,為此PacBio公司推出了采用Hifi (High fidelity reads)模式的測序方式,使得測序錯誤率大為降低,但是測序長度有所限制(15-20 kb)。為得到高質量基因組,可采用不同測序平臺測序,再進行基因組組裝,以實現優勢互補。本研究采用Nanopore長片段測序平臺的測序結果進行厚葉木蓮的基因組初步組裝,再用MGI短片段測序平臺的測序結果對組裝的基因組進行糾錯,得到了較完整和準確的組裝結果,但是還沒有組裝到染色體級別,后期本研究將繼續采用Hi-C(High-throughput chromosome conformation capture)測序方式進一步對厚葉木蓮基因組進行測序,測得的結果可用于進一步組裝以完善厚葉木蓮基因組。

在采用Nanopore測序數據進行厚葉木蓮基因組組裝過程中,本研究采用NextDenovo 2.3.1軟件進行序列拼接,這一軟件采用先糾錯再組裝(Correction then assembly)的模式進行基因組組裝,保證了組裝序列的連續性和正確性,優于其他組裝軟件。同時本研究利用多種方法去除了原始組裝序列中的冗余序列并進行后期糾錯:其中racon v1.5.0和hapo-G v1.3.2分別使用Nanopore測序獲得的長片段和MGI測序獲得的短片段進行組裝基因組內堿基和插入缺失錯誤的糾錯,polypolish v0.5.0主要針對組裝序列中的同聚體長度(Polypolish-length,如“AAAAAA”這種重復序列的長度)進行糾錯,3個糾錯軟件的使用可以很好地保證組裝基因組的準確性。

在基因組注釋環節,本研究利用EDTA和RED軟件相互組合的方法進行基因組中重復序列的查找。EDTA主要用于基因組中轉座子的查找,而RED查找包括轉座子在內的所有重復序列類型,如微衛星體(Microsatellite)等,兩者結合保證了重復序列查找的全面性,但RED分析的結果并沒有對不同重復序列進行歸類。本研究首先使用BRAKER2,結合轉錄組數據和其他物種的蛋白質序列預測厚葉木蓮基因組中的基因,得到初步基因序列注釋結果。由于這一結果還存在較多錯誤預測,本研究進一步利用funannotate軟件對初步預測的結果進行整合。Funannotate是個軟件整合工具,被編譯為一個管段(Pipeline)流程。該軟件可以進行不同方式的基因預測(包括de novo注釋、轉錄組注釋和同源蛋白質注釋),再結合其他分析軟件的結果(本研究中為BRAKER2基因預測結果),可以獲得統一、可靠的高質量基因預測結果。

3.2 厚葉木蓮的比較基因組學研究

本研究中厚葉木蓮的基因組組裝大小為2 092 298 891 bp(約2.09 Gb),大于木蘭科的日本厚樸(Magnoliahypoleuca,修正名為Houpoeaobovata)基因組(約1.64 Gb)[58]、厚樸(Magnoliaofficinalis,修正名為H.officinalis)的基因組(約1.68 Gb)[59]和鵝掌楸的基因組(約1.74 Gb)[60],小于望春玉蘭的基因組(約2.22 Gb)[29]?；凇癳mbryophyta_odb10.2020-09-10”單拷貝基因庫,厚樸基因組的組裝完整性評估為86.20%的BUSCO單拷貝基因能被完整匹配到,日本厚樸基因組的組裝完整性為98.6%,望春玉蘭基因組的組裝完整性為95.7%(上述數值為本研究利用望春玉蘭基因組分析的結果),鵝掌楸基因組的組裝完整性為98.8%(此值為本研究利用鵝掌楸基因組重新分析的結果),而厚葉木蓮的BUSCO組裝完整性為98.8%,與鵝掌楸的結果相同,但高于厚樸、日本厚樸和望春玉蘭3個物種。從重復序列比例看,厚葉木蓮基因組中重復序列占基因組的76.5%,高于望春玉蘭(66.48%)[29]、日本厚樸(64.54%)[58]以及鵝掌楸(63.81%)[60],但小于厚樸(81.44%)[57],厚葉木蓮基因組重復序列比例在已有木蘭科物種基因組重復序列比例的范圍內。上述結果說明,雖然本研究報道的厚葉木蓮基因組還是草圖狀態,沒有組裝到染色體級別,但是組裝已經很完整,結果可靠。

對基因家族的研究表明,厚葉木蓮中與次生代謝相關的基因家族顯著擴張,如萜類物質,這與日本厚樸、望春玉蘭基因家族分析結果相似[29,58]?；虮磉_分析表明,望春玉蘭花中萜類相關基因的表達高于其葉子,說明萜類物質是望春玉蘭花香的主要原因[29],厚葉木蓮花是否也有類似的結果還需結合基因組進一步研究確定。除了次生代謝物質相關的基因家族在厚葉木蓮基因組中有顯著擴張外,本研究還發現與木質部、韌皮部發育相關的基因家族,以及與光合作用[光合作用光反應(Photosynthesis light reaction,GO:0019684)、光合作用光收獲(Photosynthesis light harvesting,GO:0009765)以及KEGG中的光蛋白(Photosynthesis proteins)]、溫度[熱反應(Response to heat,GO:0009408)、溫度刺激反應(Response to temperature stimulus,GO:0009266)]等相關的基因家族也有擴張。厚葉木蓮生長在較高海拔的地區,木材致密通直[26],這些擴張的基因家族是否與其適應高海拔的溫度、光照、強風等生境有關還需進一步研究。厚葉木蓮基因組中的擴張基因家族中的一部分與肌動蛋白絲(也稱“微絲”)相關,而肌動蛋白是植物細胞骨架的主要元素,相關研究表明它與植物抗真菌功能密切相關[61,62]。

對日本厚樸、望春玉蘭的研究均表明木蘭科都共同經歷了一次近期的全基因組重復事件[29,58]。厚葉木蓮有著與日本厚樸相似的全基因組重復基因、串聯重復基因和近端重復基因比例。與全基因組重復相關的基因占厚葉木蓮基因組基因的比例為12.6%,日本厚樸為13.4%[58];串聯重復基因在厚葉木蓮基因組基因的比例為8.8%,日本厚樸為7.6%[58];近端重復基因在厚葉木蓮基因組基因的比例為8.2%,日本厚樸為9.4%[58]?；虼撝貜秃徒酥貜褪腔蛐纬傻闹匾绞?與物種適應性密切相關[57]。日本厚樸中串聯重復和近端重復基因的功能主要與苯丙烷類、萜類、類黃酮的生物合成等有關[58],這與厚葉木蓮這兩類基因的研究結果相同,但厚葉木蓮可產生更多的次生代謝合成物,如苯并惡唑嗪酮類化合物、生物堿、聚酮等。厚樸基因組中,萜類合成相關的基因也呈串聯重復狀況[59]。值得注意的是日本厚樸原產日本,具有很好的抗寒性,研究發現苯丙烷類生物合成與日本厚樸抗寒性相關。

對厚葉木蓮中與全基因組重復相關基因的富集分析表明,這些基因與植物最基本的生長發育密切相關,如植物晝夜節律相關基因[63]、與碳代謝相關的檸檬酸循環基因[64]、參與細胞的多種生命活動(如細胞通訊、核糖體與內質網的結合、小泡運輸、蛋白質合成等)的GTP結合蛋白基因[65,66]等,說明全基因組重復在木蘭科植物適應性方面具有重要作用。

4 結論

本研究首次報道了木蘭科木蓮屬物種的基因組,這一結果對全面深入了解木蘭科物種以及厚葉木蓮的進化、適應具有重要作用。本研究組裝的厚葉木蓮基因組大小為2 092 298 891 bp,基因組序列中76.5%為重復序列。通過基因預測,在組裝的厚葉木蓮基因組中共注釋到37 900個基因,它們編碼了41 675種蛋白質。對厚葉木蓮基因組研究發現,全基因組重復與木蘭科植物進化適應性密切相關,很多與植物基本生長發育相關的基因通過全基因組重復在木蘭科得到加強;厚葉木蓮富含有與次生代謝相關的多種基因,由此產生的次生代謝物質有利于厚葉木蓮在高海拔生長以及抵御病蟲害,但其(種子)香味也使得其易受啃食傷害。因此,厚葉木蓮基因組的組裝為從遺傳角度深入了解其瀕危機制提供了可能,也為科學合理利用厚葉木蓮以及提取其次生物質作為生物醫藥提供了參考。