20230402_基于染色体水平基因组组装揭示印楝的萜类化合物的生物合成

印楝（Azadirachta indica，又称苦楝）是楝科的常绿树，印楝素是印楝属植物根、茎、叶等部位所含有的一类具有重要药理活性的化合物，具有抗肿瘤、抗病毒、抗炎等多种药理作用。由于缺乏染色体级别的组装，我们对于印楝树基因组结构和A. indica的比较基因组分析的深入理解受到了阻碍。本文采用Illumina、PacBio和Hi-C技术结合构建了A. indica高质量基因组组装，这是A. indica的第一个染色体级别的基因组组装。根据我们的组装长度，A. indica基因组大小估计为281 Mb，锚定在14条染色体上（contig N50 = 6 Mb，scaffold N50 = 19 Mb）。该基因组组装包含了115 Mb的重复元件和25,767个编码蛋白基因。我们还鉴定了次生代谢产物合成基因簇、TPS基因和CYP基因。比较基因组分析显示，大多数A. indica特异性的TPS基因和CYP基因位于第13号染色体上的萜类相关群中，这表明第13号染色体在A. indica特异性萜类合成中发挥着重要作用。我们所提供的A. indica基因组数据和基因组分析将有助于我们更好地了解A. indica中萜类的合成，并促进楝科家族的比较基因组研究。同时，印楝染色体组装可以促进对印楝素生物合成途径的深入理解，为印楝素的合成机制和药理作用的研究提供有力的支撑，具有重要的科学价值和应用前景。

A. indica的基因组组装的基因组特征。（A）A. indica的Hi-C序列在相同以及不同染色体之间的交互频率；（B）A. indica基因组组装和注释结果；（A）不同的染色体。（B）基因在染色体不同位置上的分布密度。以500 kb长度为单位，计算每个单位内的基因密度；（C-G）A. indica所有基因在不同染色体上的表达量（从外到内轨迹：茎、根、叶、果实和花）；（H，I）染色体不同位置上的重复序列分布密度以及GC碱基分布密度（J）次生代谢相关基因在染色体上分布的位置。（K）不同染色体以及相同染色体内的共线性模块。

• 次生代谢物分析结果。 (A) A. indica 和其他10个物种的次生代谢物分析。(B) 13号染色体上鉴定到的萜类相关基因的基因簇(C) 13号染色体上基因的KEGG通路富集分析的点图。(D) 13号染色体上基因的GO富集分析 (P < 0.05)。
  
• 植物基因组中编码某些特殊代谢途径的基因通常会聚集在一起。我们利用PlantiSMASH分析流程鉴定了A. indica基因组中的参与特殊代谢途径的聚集的基因聚集群。根据分析结果，在A. indica基因组中鉴定出了50个聚集群，包括692个基因。鉴定到的聚集群大小范围为27.2至1634.4 kb。此外，在其他11个物种中也检测到了41个（C. sinensis），51个（A. yangbiense），48个（T. cocoa），47个（G. raimondii），45个（A. thaliana），35个（F. vesca），33个（P. persica），30个（C. sativus），46个（V. vinifera），47个（S. lycopersicum）和29个（B. distachyon）群集（图4A）。如预期的那样，在A. indica基因组中识别到了更多与萜类物质相关的聚集群。
  
• Azadirachtin是印楝树中的三萜类化合物，对广泛的昆虫物种具有有效的杀虫活性，但对哺乳动物的毒性非常低。据报道，单萜合酶（TPS）、细胞色素P450酶（CYP450）、醇脱氢酶（ADH）、酰转移酶（ACT）和酯酶（EST）参与了印楝内酯的生物合成。在本研究中，大量编码CYP 450（78）、TPS（58）和ACT（34）的基因在次生代谢物生物合成基因簇中被鉴定出来。编码ADH和EST的基因可能分散在基因组中。萜类相关的簇主要分布在1、2、3、5、6、7、10、11、12和13号染色体上。四个萜类相关的簇（簇18-21）覆盖了约1.4 Mb的13号染色体上（图4B）。在13号染色体上的83个聚类的萜类相关基因中，有12个基因包含在A. indica特异扩增基因家族中。这些基因被认为是参与A. indica特有的萜类生物合成的潜在基因。
• KEGG结果显示13号染色体上的基因主要涉及“内质网中的蛋白质加工”、“倍半萜和三萜生物合成”以及“卵巢类固醇合成”（图4C）。此外，我们使用13号染色体上的基因进行GO富集分析时，与“倍半萜合酶活性”（GO：0010333）相关的基因呈现出较低的P值，说明“倍半萜合酶活性”是染色体13上富集的最主要的功能类别（图4D）。
  

• TPS 基因家族的特征是两个大结构域：PF01397（萜烯合酶，N 末端结构域）和 PF03936（萜烯合酶家族，金属结合结构域）。 为了研究 TPS 基因家族的特征和进化，我们在 A. indica 和其他 13 种植物基因组中鉴定了总共 512 个推定的 TPS 基因。 在 A. indica 中鉴定出 70 个推定的 TPS 基因； 这些包括 44 个包含 PF01397 和 PF03936 结构域的 AziTPS 基因、9 个包含 PF01397 结构域的 AziTPS 基因和 17 个包含 PF03936 结构域的 AziTPS 基因。 与其他植物相比，A. indica (N = 70) 包含最多的 TPS 拷贝，其次是 C. sinensis (N = 49) 和 A. yangbiense (N = 57)（图 5A)
• 使用来自A. indica和其他13种植物的403个TPS基因进行系统发育分析，其中有6个已知的TPS基因属于TPS-a、b、c、e、f和g亚家族，其余的TPS基因长度太短，无法进行有意义的比对。如图5B所示，6个亚家族的拓扑结构与先前的研究结果相似。在40个用于系统发育分析的AziTPS中，有15个、13个、4个、3个、1个和4个AziTPS基因分别属于TPS-a、b、c、e、f和g亚家族。在A. indica中，TPS-a和-b亚家族是主要亚家族，在系统发育分析中占总AziTPS基因数量的约37.5%和32.5%。此外，该研究还使用系统发育分析和55%相似性的阈值，鉴定出了可能是A. indica特有的TPS基因，并发现了9个A. indica特有的AziTPS基因。有趣的是，其中7个特定的AziTPS基因（Indica_007028、Indica_007047、Indica_007053、Indica_007068、Indica_007070、Indica_007072和Indica_007143）位于第13号染色体的萜烯相关聚类群（聚类群18、19和20）中。
• 我们进一步研究了 A. indica 中 TPS 基因的表达模式。 来自 A. indica 五个组织的转录组数据集是从我们之前的研究（Wang et al., 2020）中获得的，并在本研究中重新映射到染色体水平的基因组组装。 超过 88% 的 RNAseq 读数被唯一映射到所有样本的基因组组装。 在A. indica 的五个组织中共检测到 27 个 TPS 基因表达的转录本。 大多数检测到的转录物都表现出空间特异性表达模式（图 5C）。 有九个、一个、四个、三个和四个基因分别在花、果实、根、叶和茎中专门表达。 七个基因（AziTPS30、-48、-5、-57、-26、-63 和-50）主要在一个或两个组织中表达。

• 从所有 A. indica 和其它13个植物基因组中共鉴定出 3,657 个 CYP 基因（图 6A ）。 A. indica 基因组中共有 355 个 CYP 基因。 此外，用 A. indica 的 157 个全长 CYP（450 < 长度 < 600）蛋白质序列构建系统发育树结果如图 6B 所示，系统发育树分为两个主要进化枝：A 型（49%；77/157）和非 A 型（51%；80/157）； 并进一步聚集成九个亚型。 71 氏族是最大的氏族，由 49% (77/157) 的成员组成； 72亚型、74亚型、85亚型、86亚型分别包含18、4、28、25个CYP基因。剩下的51亚型、710亚型、711亚型、727亚型均只包含一个CYP基因。
  
• 为了鉴定A. indica 特有的 CYP 基因，作者利用氨基酸序列比对构建了一个包含2807个（长度在450到600之间）来自A. indica和其他13个植物基因组的CYP基因的系统发育树，并使用55%相似度进行过滤。最终鉴定出了6个A. indica 特有的CYP基因。与TPS基因类似，这5个CYP基因（Indica_007272、Indica_007273、Indica_007276、Indica_007277和Indica_007278）中的五个位于第13号染色体的萜类次生代谢产物生物合成基因簇。这些位于第13号染色体上萜类次生代谢产物生物合成基因簇的特异性TPS和CYP可能参与印楝素的特异性生物合成。
  
• 作者还研究了A. indica CYP基因在不同组织（果实、花、根、茎和叶）中的表达模式。221个CYP基因的转录本在不同的组织中被检测到（图6C）。果实、茎和叶中有更多的高表达CYP基因，而花和根中则较少。果实、茎和叶中的高表达CYP基因分别为83、88和97个。在含有高azadirachtin A含量的组织（果实和叶）中高表达的CYPs更可能参与azadirachtin的生物合成。此外，A. indica特异性的AziCYP256（Indica_007272）和AziCYP8（Indica_007273）在果实和花中高表达。
  

===================================

DNA制备和基因组测序方法

作者从海南大学留芳园采集新鲜的印楝树叶子作为基础材料，并利用DNAsecure植物提取试剂盒从中提取出基因组DNA。研究中采用了Illumina HiSeq X Ten平台进行双端测序和PacBio RSII平台进行长读测序。其中Illumina双端文库插入片段大小为270 bp，PacBio文库插入片段大小为20 kb。

基因组组装方法

作者使用Canu v2.0软件对原始的PacBio测序reads进行矫正和组装，得到了886个contig，N50大小约为6 M。然后使用RACON软件对组装的contig进行一轮修正，并利用Illumina短读对polished contigs进行两轮修正。最终，870个contig通过RACON和Pilon软件进行了错误校正。A. indica基因组大小是通过流式细胞术法估计的。

Hi-C辅助染色体水平组装

为了进行Hi-C测序，我们生成了一个150bp的配对末端文库，并在Illumina HiSeq X Ten平台上进行了测序。使用默认参数的Bowtie2软件将干净的reads映射到A. indica基因组上。使用HiC-Pro v2.11.1软件将Hi-C测序reads映射到组装的草稿基因组上并检测有效的接触点。然后，我们使用ALLHIC v0.9.12软件根据有效接触点之间的关系将contig聚类成染色体级别的scaffold。

基因组完整性的评估

首先使用基于OrthoDB1数据库的BUSCO软件评估了基因组装配的准确性和完整性。随后使用LTR_Retriever计算了染色体的LAI分数。接下来，将来自A.indica 的5个组织的转录组NGS序列用Hisat2比对到染色体上并统计最终比对率。同时，也使用Bowtie2将基因组NGS短读序列映射到了染色体上并统计最终的比对率。最后，使用Minimap22和dotPlotly.3对我们的装配的染色体与已发表的contig水平基因组进行了线性分析比较。

重复序列的鉴定

作者使用RepeatModeler v1.0.8 (Price et al., 2005)和RepeatMasker v4.0.7 (Tarailo-Graovac和Chen，2009)两种软件鉴定重复序列。使用LTRharvest 和LTR_Finder 鉴定了A. indica的LTRs。随后使用LTR_retriever 整合了LTRharvest和LTR_Finder的结果。接下来使用RepeatModeler来鉴定interspersed repeats，并将LTR-retriever鉴定的重复序列与RepeatModeler鉴定的重复序列整合成一个本地重复序列数据库。最后作者使用RepeatMasker利用本地构建好的重复序列数据库和Repbase数据库基于同源性比对的方法鉴定出A. indica基因组中的重复序列。

非编码RNA的鉴定

通过使用tRNAscan-SE 对可靠的tRNA进行搜索。使用INFERNAL 在Rfam 数据库中搜索小核RNA (snRNA)和微小RNA (miRNA)。

基因预测

使用三个代表性物种（包括Citrus sinensis，Theobroma cacao和Acer yangbiense）的基因组进行同源性注释。使用TBLASTN软件将这些物种的蛋白质序列与A. indica基因组序列进行比对，E值≤1e-5。根据TBLASTN结果，使用GeneWise预测确切的基因结构。使用Cufflinks 根据RNA-seq数据初步识别基因结构。使用Augustus 和SNAP对染色体进行从头预测。使用MAKER软件对上述三种方法预测得到的所有基因进行汇总整合。

基因注释

将预测的基因序列与四个蛋白质数据库，包括NR、InterPro、Swiss-Prot和EggNOG进行比对和注释。使用PlantiSMASH从A. indica基因中预测和注释与代谢相关的基因群。使用DRAGO预测A. indica基因中的疾病抗性基因。

鉴定萜烯合成酶和细胞色素P450家族相关基因并构建系统发育树

使用HMMER 对鉴定到的基因与Pfam-A数据库进行比对，E值设定为1e-5。PF01397（萜烯合酶，N端结构域）和PF03936（萜烯合酶家族，金属结合结构域）域被用来识别TPS基因家族成员。用于进行系统发育分析的数据集包括来自A. indica和其他13种植物的403个潜在TPS，以及属于TPS-a（AAX16121.1）、b（AAQ16588.1）、c（AAD04292.1）、e（Q39548.1）、f（Q93YV0.1）和g（ADD81294.1）亚家族的6个已报道的TPS基因。PF00067（细胞色素P450）被用于识别CYP基因家族成员。通过氨基酸长度筛选出潜在的CYP基因（450<长度<600），用于进行系统发育分析。系统发育树主要有MEGAX软件基于TPS和CYP蛋白质序列的比对结果而构建的。