毛茛科黑种草Nigella damascena巨型基因组--nature communications

The Nigella damascena genome provides insights into transposable element-driven genomic gigantism and trait evolution

黑种草基因组解析揭示转座子驱动的基因组巨型化与性状演化机制

摘要

被子植物中已多次演化出以转座子（TE）为主体的巨型基因组。然而，转座子在巨型基因组演化中发挥的作用仍尚不明确。本研究通过祖先基因组大小重建、全基因组测序，以及比较基因组与转录组分析，揭示毛茛科（Ranunculaceae）基因组巨型化的发生过程、驱动因素及其演化后果。研究发现，大马士革黑种草（Nigella damascena）的二倍体巨型基因组（10.84 Gb，约为耧斗菜基因组的 35 倍）起源于约 1.46 Gb 的祖先基因组，在 6000 万年的漫长时间里通过转座子持续积累形成。基因区域内的转座子插入产生了内含子超长的基因、编码序列发生改变的基因，以及黑种草属特有的转座子起源基因；这类基因合计约占大马士革黑种草蛋白编码基因的 20%。转座子介导的调控改变、基因复制与丢失，很可能是大马士革黑种草精巧花瓣、合生心皮以及特化次生代谢产物演化的分子基础。上述结果阐明了转座子如何驱动基因组巨型化并塑造物种特化性状，加深了我们对基因组大小演化的认识。

引言

被子植物的基因组大小差异极其悬殊，变异幅度可达 2400 倍以上：记录的最小基因组为狸藻科的Genlisea tuberosa（61 Mb），最大基因组为黑药花科的日本重楼（Paris japonica，148.8 Gb）。基因组大于 10 Gb 的巨型基因组物种，已见于被子植物 24 个目下至少 43 个科，说明巨型基因组在被子植物中经历过多次独立演化。普遍认为，基因组巨型化深刻影响对应类群的生存，还有可能促成植物多样形态与功能的起源。

理论上，多倍化（即全基因组复制 WGD）与重复序列扩增均可驱动基因组扩张。但越来越多的基因组证据表明，巨型基因组的主体通常是重复序列，尤以转座子（TE）占比最高，多数物种中转座子占基因组 80% 以上，例如石蒜科的洋葱、菊科的菊芋、芍药科的滇牡丹、百合科的川百合、黑药花科的云南重楼。该类反复出现的特征说明，重复序列扩增在巨型基因组形成过程中发挥重要作用。

大量案例研究证实，转座子插入基因组特定区域可以调控基因表达与功能、引发基因丢失或复制，甚至催生新基因，进而推动植物性状演化。例如，玉米驯化过程中，反转录转座子 Hopscotch 插入TEOSINTE BRANCHED1（TB1）基因上游，造成该基因过表达，使玉米分蘖数显著减少；反转录转座子 Dasheng 插入Wx等位基因第 9 外显子，产生功能丧失突变，形成糯稻品种；反转录转座子 Rider 介导新基因SUN诞生，是番茄果实由圆形演变为长形的分子基础。

已有研究报道，在大蒜、滇牡丹、王百合等拥有巨型基因组的物种中，相当比例的蛋白编码基因内部含有转座子序列。在滇牡丹中，转座子影响油脂生物合成、花性状（花色、花香）相关基因的演化。尽管已有上述进展，但在漫长演化历程中转座子如何塑造巨型基因组，转座子扩增在多大程度上影响基因演化并最终塑造特化性状，仍知之甚少。

毛茛科是研究上述科学问题的理想研究体系。一方面，毛茛科是被子植物中基因组大小跨度极大的类群。该科部分属，如黑种草属、毛茛属、獐耳细辛属，拥有大甚至巨型基因组，染色体长且弯曲，称为毛茛型（R-型）染色体。这和近缘类群形成鲜明对比：耧斗菜属为小型肾形的唐松草型（T-型）染色体；黄连属为细小纤长的黄连型（C型）染色体。毛茛科基因组分化显著，系统发育关系解析完善，已发表参考基因组的基因组大小跨度约 300 Mb～6.0 Gb；随着巨型基因组数据的补充，毛茛科成为解析基因组巨型化演化的关键研究材料。

大小不一的彩色圆圈代表分支节点处推断的祖先基因组大小，以及末端现存物种的基因组大小。右侧代表性物种的染色体图改编自 Gregory 等文献，空心圆圈代表染色体末端的卫星序列。绿色星号标记毛茛目与翠雀族起源之前发生的全基因组复制（WGD）事件；结合后续共线性分析、共线性旁系同源与直系同源基因的 Ks 分布比较，估算了这些 WGD 事件的发生时间。Mya：百万年前。原始数据见源数据文件。

另一方面，毛茛科物种演化出了独特的花部性状与植物化学特征。例如，花的对称性可为辐射对称或两侧对称；花器官排列方式有轮状排列与螺旋状排列；花瓣与退化雄蕊可存在或缺失；花瓣形态可呈扁平、具距或盾状；心皮可离生或合生。毛茛科物种还以丰富、类群特异性的植物次生代谢产物闻名，包括原小檗碱型生物碱、二萜生物碱、毛茛苷以及 β- 榄香烯。研究人员以耧斗菜属、唐松草属、翠雀属和黑种草属为模式材料，对花发育与花瓣形态建成相关的遗传调控通路开展了大量研究。已有数据表明，该科内花器官身份基因的功能整体保守；花瓣的特化与多样化，源于背腹轴基因的表达重排，以及原本不在花瓣中表达的基因被共选择（共 opt）进入花瓣发育通路。尽管已有上述机制层面的认识，但毛茛科特化性状演化的基因组基础仍不清楚。

本研究基于可靠的时间标定系统发育框架，重建毛茛科祖先基因组大小，并追溯巨型基因组的演化轨迹。我们组装了大马士革黑种草（一种典型二倍体巨型基因组物种）的染色体水平基因组，结合多个基因组大小差异显著的毛茛科物种（蓝耧斗菜 310 Mb、黄连 0.94 Gb、草甸毛茛 2.69 Gb、云南乌头 5.76 Gb）开展大规模比较基因组与转录组分析，揭示转座子扩张是驱动基因组膨胀的核心因素，及其对基因演化的广泛影响。此外，我们解析了大马士革黑种草多个特化性状演化的基因组基础，包括二唇形花瓣特化、心皮后生融合、β- 榄香烯合成以及原小檗碱生物碱合成通路丢失。本研究结果不仅加深了对基因组大小演化的理解，也阐明转座子在塑造巨型基因组、塑造物种特化性状中的关键作用。

结果

毛茛科基因组大小的演化

为解析毛茛科基因组大小的演化过程，本研究整合前人研究构建时间标定系统发育树，基于该框架开展祖先基因组大小重建（图 1，附图 1）。数据集一共包含 68 个物种：50 个毛茛科物种，18 个其他毛茛目类群作为外类群；基因组大小数据来自流式细胞仪测定或已发表文献（附图 1、补充数据 4，参考文献见其中）。以系统发育树与基因组大小为输入，采用 phytools 包中anc.ml函数，基于布朗运动（BM）模型重建祖先状态（补充数据 4）。结果显示：毛茛科最近共同祖先基因组大小为 1.46 Gb；核心毛茛科（排除 3 个基部族）的最近共同祖先基因组大小为 2.61 Gb，二者均与现存具有 C 型染色体物种的基因组规模接近（图 1）。核心毛茛科分化之后，基因组大小朝两个截然相反方向演化：具有 T 型染色体的人字果族，其最近共同祖先基因组剧烈收缩至 0.61 Gb；其余具有 R 型染色体类群的最近共同祖先基因组扩张至 3.56 Gb（图 1）。该 3.56 Gb 的祖先基因组随后继续扩张，在不同属中独立演化出大型甚至巨型基因组，例如獐耳细辛属、毛茛属、菟葵属、Anemonopsis、黑种草属、侧金盏花属（图 1）。以黑种草属演化支为例：基因组由祖先状态 3.56 Gb 扩张至该属最近共同祖先的 6.83 Gb，多数现存物种基因组持续扩张，大马士革黑种草与东方黑种草基因组甚至超过 10 Gb（图 1）。以上结果说明，毛茛科演化过程中基因组大小经历了剧烈扩张或收缩；现存类群的巨型基因组是多次独立、逐步扩张形成的。

大马士革黑种草的巨型基因组及其核心特征

为探究毛茛科基因组扩张的驱动因素与演化后果，本研究对二倍体巨型基因组代表物种大马士革黑种草开展基因组测序与特征解析。利用 274.82 Gb HiFi 长读长，得到初始组装版本，基因组大小 10.84 Gb，包含 3192 条 contig，contig N50=11.58 Mb；该组装大小与流式细胞仪（10.32 Gb）、k-mer 估算（10.82 Gb）结果基本吻合（附图 1、补充数据 4）。再利用 247.04 Gb 有效 Hi-C 交互数据，将 contig 锚定至 6 条假染色体，挂载序列占组装基因组的 95.86%（10.39 Gb）（图 2，附图 1、补充数据 4）。该基因组组装完整性与连续性优异：Illumina 回比对率高于 99.25%，BUSCO 完整度 98.50%，LTR 组装指数 LAI 约为 17（附图 1、补充数据 4）。

a 大马士革黑种草基因组特征。a 轨道代表组装得到的 6 条染色体；每条染色体上红色方框为组装的着丝粒区域，黑色方框为端粒区域。b 轨道代表 GC 含量。c-g 轨道分别为蛋白编码基因、非编码 RNA、转座元件、Gypsy、Copia 的密度。bg 轨道基因组特征均采用30 Mb 滑动窗口、步长 3 Mb计算。内部连线代表共线性区块。

b 7 个毛茛目物种的全基因组共线性。彩色连线连接共线性区块，每种颜色对应大马士革黑种草的一条染色体；一个代表性共线性区块以绿色高亮标注。带编号白色椭圆代表染色体。

c 云南乌头、大马士革黑种草、虞美人的微共线性比较。蓝色矩形代表正链蛋白编码基因，绿色矩形代表负链蛋白编码基因；灰色连线连接共线性基因对，其中一组分别以紫色、棕色高亮显示。

d 经演化校正后的共线性旁系同源与直向同源基因 Ks 分布。红色箭头指示大马士革黑种草内部共线性旁系同源基因的 Ks 峰；灰色箭头标记大马士革黑种草与虞美人之间共线性直向同源基因的 Ks 峰。 物种缩写：Nida，大马士革黑种草；Acvi，云南乌头；Raca，草甸毛茛；Aqco，蓝耧斗菜；Coch，黄连；Eppu，Eranthis pubescens；Parh，虞美人。

本基因组共注释得到35 747 个蛋白编码基因，重复序列总长 9.61 Gb（表 1，补充数据 4）。进一步与云南乌头、草甸毛茛、蓝耧斗菜、黄连的基因组对比发现：大马士革黑种草多项核心基因组指标与其余 4 个物种大体相当，甚至数值更小，包括蛋白编码基因数量、平均编码序列（CDS）长度、平均外显子长度、单基因平均外显子数目（表 1）。 与之相反，大马士革黑种草的平均基因长度（9715 bp）、平均内含子长度（3016 bp）、重复序列含量（9.61 Gb，占基因组 88.60%）显著更高（表 1）。具体来看，大马士革黑种草重复序列中转座元件总长 9.59 Gb，分别是其余 4 个物种转座元件总量的 1.94 倍、4.14 倍、69.76 倍、14.49 倍（表 1，补充数据 4）。

我们随后对黑种草（Nigella damascena）开展基因组内共线性分析，鉴定得到 96 个基因组内共线性区块，证实该物种曾发生一次明显的全基因组复制（WGD）事件（图Fig2，补充图MOESM1）。后续的种间共线性分析结果显示：黑种草与滇川翠雀（Aconitum vilmorinianum）之间存在一对二的共线性区块；而黑种草与另外 3 种毛茛科物种 —— 肾叶唐松草（Thalictrum cassubicifolius）、蓝花耧斗菜（Aquilegia coerulea）、中国金莲花（Caltha chinensis），以及小檗科的淫羊藿（Epimedium pubescens）、罂粟科的虞美人（Papaver rhoeas）之间，均为一对一共线性区块。这说明毛茛目类群共享一次古老的全基因组复制事件，而滇川翠雀额外发生过一次更近的全基因组复制事件（图Fig1、Fig2，补充图MOESM1）。

同义替换率（\\(K_S\\)）分析表明，每个物种基因组内的共线性旁系同源基因均呈现单一峰值，且\\(K_S\\)＞1.0（补充图MOESM1），进一步佐证毛茛目物种共享这一古老全基因组复制事件。经进化速率校正后，两次全基因组复制事件的发生时间分别约为 1.2812 亿年前（校正后\\(K_S\\)峰值 = 1.08）、5060 万年前（校正后\\(K_S\\)峰值 = 0.43），两次事件均早于毛茛目以及翠雀族的起源时间（图Fig1、Fig2）。综上可知，黑种草的基因组扩张源于大量转座元件（TE）的积累，而非全基因组复制；本研究检测到的这次 WGD 事件，发生在毛茛目分化之前。

黑种草演化历程中转座元件扩张的模式与动态特征

为阐明转座元件如何塑造巨型基因组，本研究对比了黑种草、滇川翠雀、肾叶唐松草、蓝花耧斗菜、中国金莲花的转座元件积累模式。在基因区范围内，黑种草的转座元件总长 162.45 Mb，与滇川翠雀（132.59 Mb）大体相近，但显著高于肾叶唐松草（33.01 Mb）、蓝花耧斗菜（11.45 Mb）、中国金莲花（36.57 Mb）（图Fig3，补充数据MOESM4）。在基因间区，黑种草共积累 9.40 Gb 的转座元件，分别是滇川翠雀、肾叶唐松草、蓝花耧斗菜、中国金莲花的 1.97 倍、4.16 倍、80.44 倍、15.78 倍（图Fig3，补充数据MOESM4）。

五大主要转座元件类别中，黑种草的长末端重复反转录转座子（LTR，9.18 Gb）与 DNA 转座子（387.30 Mb）的数量、序列长度都远高于其余 4 个物种；相较于滇川翠雀、肾叶唐松草、蓝花耧斗菜、中国金莲花，LTR 含量依次为 2.07、4.33、112.27、17.66 倍，DNA 转座子依次为 2.32、10.61、35.56、8.19 倍（图Fig3，补充数据MOESM4）。以上结果说明：基因间区转座元件的大规模扩增，叠加基因区内转座元件的适度积累，是黑种草基因组巨型化的主要驱动因素，其中 LTR 反转录转座子和 DNA 转座子贡献最高。

a 5 个基因组的总转座元件（TE）含量。b 代表基因组区域，展示 5 个物种的转座元件密度与类型；局部放大两个子区域，凸显黑种草基因区与基因间区的转座元件积累情况。灰色连线代表共线性基因对。c 5 个基因组中 5 类转座元件的序列长度。di 5 个物种全部转座元件（d）或特定类型转座元件（ei）的估算插入时间。每一列代表某一插入时间下，全部转座元件（d）或特定类型转座元件（ei）占 5 个物种转座元件总含量的比例；各物种比例分别独立计算。Mya，百万年前。

我们进一步比较了 5 种毛茛科物种的转座元件插入动态与清除效率。本研究采用两种策略估算转座元件的插入时间与积累动态。仅对完整 LTR 反转录转座子开展分析，结果显示 5 个物种在过去 1000 万年内均发生过转座元件近期爆发事件（补充图MOESM1）。与之不同，当纳入全部已注释转座元件进行分析时，具有不同染色体类型的物种呈现出截然不同的积累模式。 拥有 R型染色体的物种（黑种草、滇川翠雀、肾叶唐松草）自 6000 万年前起转座元件持续积累，该时间大致与该类染色体类群的演化辐射时间重合（图Fig1、Fig3）。具有 C型染色体的中国金莲花，转座元件积累发生时间相对更晚；而拥有 T型染色体的蓝花耧斗菜，并未检测到显著的转座元件爆发事件（图Fig3）。 在黑种草中，约 92% 的转座元件插入发生在黑种草属与其姐妹类群（翠雀族，约 5900 万年前）分化之后；约 2700 万年前的转座元件爆发与黑种草属早期分化时间相吻合，反映出转座元件以类群特异性积累为主（图Fig1、Fig3）。此外，不同类型转座元件对该物种基因组巨型化的贡献存在差异：LTR 反转录转座子与短散在核元件（SINE）的积累模式与全部转座元件整体趋势基本一致；而 DNA 转座子、未分类转座元件、长散在核元件（LINE）的爆发峰值出现更早，分别为 3700 万年前、4400 万年前、5100 万年前（图Fig3）。

通过计算完整 LTR 反转录转座子与单拷贝 LTR + 截短型 LTR总和的比值，评估转座元件清除效率。单拷贝 LTR 与截短型 LTR 是两种主要转座元件清除机制的产物，即不等同源重组与双链断裂（DSB）介导的序列缺失。结果显示，除肾叶唐松草（比值约 1:8）以外，其余 4 个物种均表现出极高的转座元件清除效率：黑种草约 1:195、滇川翠雀约 1:199、蓝花耧斗菜约 1:70、中国金莲花约 1:117（补充数据4），该特征与基因组大小无关。综上表明：黑种草的基因组巨型化由长期、渐进式的转座元件积累驱动；而蓝花耧斗菜基因组收缩，则是转座元件活性低与高效 DNA 清除共同作用的结果。

转座元件扩增对基因演化的整体影响

为评估大规模转座元件扩增对基因演化造成的影响，本研究在黑种草、滇川翠雀、肾叶唐松草、蓝花耧斗菜、中国金莲花基因组中筛选携带转座元件的基因。黑种草基因组中共鉴定到 32863 个携带转座元件的基因（占全部基因的 91.93%）；其中 CDS 区含有转座元件的基因 11339 个，内含子区含有转座元件的基因 12795 个，基因侧翼区（翻译起始位点上游 1 kb 至翻译终止位点下游 1 kb）含有转座元件的基因 31511 个（图Fig4）。 不同物种之间，各区域内含转座元件的基因占比、CDS 区与侧翼区中转座元件的中位数数量和长度大体相近（图Fig4）。但黑种草基因的内含子区含有数量更多、长度更长的转座元件，每个基因内含子中转座元件数量中位数为 12 个，长度中位数为 3658.0 bp（图Fig4）。黑种草共有 4194 个基因拥有超长内含子（内含子区内转座元件总长度≥10 kb），远高于肾叶唐松草（624 个）、蓝花耧斗菜（26 个）与中国金莲花（295 个）（补充数据4）。

a 5 种毛茛科物种中，CDS 区、内含子区、侧翼区（翻译起始位点上游 1 kb、翻译终止位点下游 1 kb）携带转座元件的基因数量。括号内数值代表含转座元件基因占该物种全部蛋白编码基因的比例。 b、c 5 个物种基因的 CDS 区、内含子区、侧翼区内转座元件数量中位数（b）与长度中位数（c）。在 b、c 中，箱体代表第一、三分位数，横线为中位数；须延伸至最小值与最大值，离群点已剔除。b、c 中各组基因数量与 a 图对应。 d 超长内含子对基因表达的影响。 e 两个基因示例（E2F1：黑种草上调基因；UGE5：黑种草下调基因），展示这 3 个物种中的基因结构、转座元件分布与表达模式。数据为 3 个生物学重复的平均值 ± 标准差。 f 转座元件插入对编码序列的影响。 g 示例展示转座元件介导的外显子片段插入（KANADI 2，KAN2）以及终止密码子改变（Nida5G03533），造成黑种草编码序列与另外两个物种发生分化。红色圆圈标记引起蛋白序列改变的转座元件插入位点。 h 黑种草中两个黑种草属特有、起源于转座元件的基因示例（Nida3G05031、Nida4G00472）。所有示例中基因结构按 5’→3’方向展示，箭头代表转录方向。 Nida：黑种草；Acvi：滇川翠雀；Raca：肾叶唐松草；Aqco：蓝花耧斗菜；Coch：中国金莲花。源数据已作为源数据文件提供。

为解析内含子区转座元件积累对基因表达的影响，本研究开展比较转录组分析。结果发现，在黑种草中，这 4194 个带有超长内含子的基因，在检测的 10 种器官（根、茎、叶、花、果实、种子、萼片、花瓣、雄蕊、心皮）中，表达水平显著高于总长度更短的基因（总长度＜10 kb），说明富集转座元件的长内含子通常与转录活性提升相关。 本研究选取黑种草中 2538 个低拷贝（1≤拷贝数≤2）、带有超长内含子的基因，对比其在蓝花耧斗菜和中国金莲花中内含子转座元件较短（内含子转座元件总长≤1 kb）的同源基因：其中 46.38%（1177 个）基因表达模式保守，其余基因在黑种草与另外两个物种之间表达分化。分化基因包含 627 个上调基因、648 个下调基因，还有 86 个基因呈现器官依赖性表达变化。 作为对照，比较黑种草中 3609 个内含子转座元件较短的基因及其在蓝花耧斗菜、中国金莲花中的同源基因，仅 24.05%（868 个）表达模式保守。卡方检验显示，富集转座元件的基因，其表达保守性显著高于对照组（P＜0.01）。上述结果表明，转座元件扩张造成内含子增大，而内含子增大与转录活性增强、基因表达保守性升高密切相关。

为探究编码区转座元件插入带来的演化效应，本研究在毛茛科系统发育框架下，鉴定黑种草 11339 个含有转座元件来源编码序列基因的同源基因。序列比对结果显示，其中 415 个基因在其他毛茛科物种中存在同源基因，但编码的蛋白多出一段额外氨基酸序列：331 个基因为黑种草特有，84 个为黑种草属共有。 这类编码序列改变主要由起始密码子变化（116 个基因）、终止密码子变化（179 个基因），或两者同时改变（6 个基因）引起；此外外显子片段插入（57 个基因）、外显子化（21 个基因）、内含子保留（36 个基因）也有贡献。GO 富集分析显示，这些基因主要参与解剖结构发育、叶绿素生物合成过程调控。 本研究还鉴定得到 2473 个黑种草属特有基因，其编码序列大部分或全部起源于转座元件；其中 1548 个仅存在于黑种草，925 个仅在黑种草属物种间共有。这两组基因中，分别有 708 个（45.74%）和 393 个（42.49%）在至少一个检测器官中表达；708 个里面的 245 个、393 个里面的 137 个表现出器官特异性或偏好性表达。综上，编码区转座元件通过重塑原有基因、产生黑种草属新基因，推动了黑种草的基因演化。

转座元件对新形态性状演化的贡献

黑种草花的典型特征是结构高度复杂的二唇形花瓣以及合生心皮雌蕊。为解析该性状演化的基因组基础，以及转座元件在其中的潜在作用，本研究首先对比了黑种草、蓝花耧斗菜、中国金莲花的花发育 MADS-box 基因。 基因表达谱分析结果显示，大部分同源基因，包括花瓣身份基因APETALA3-3（AP3-3）、心皮身份基因AGAMOUS2（AG2），在不同物种间具有保守的时空表达模式。 黑种草中有 4 个基因表达模式与另外两个物种不同：NidaAP3-1仅在花瓣 / 雄蕊表达，NidaAP3-2仅在雄蕊表达；NidaAG1和NidaSEPALLATA2（NidaSEP2）表达量最低。进一步序列分析表明，这 4 个基因的蛋白序列与外显子 – 内含子结构保守，但NidaAP3-1、NidaAG1、NidaSEP2的启动子与内含子区域，转座元件密度显著高于它们的同源基因。

a 黑种草、蓝花耧斗菜、中国金莲花 4 个花发育阶段的花 MADS-box 基因表达谱。左侧展示 3 个物种成熟花及其花图式。天蓝色：萼片（Se）；紫色：花瓣（Pe）；橙色：雄蕊（St）；粉色：退化雄蕊；深绿色：心皮（Ca）。 b、c 二唇形花瓣复杂化的潜在机制。b NidaAP33、AqcoAP33、CochAP33 共有靶基因与特有靶基因的数量，每组列出代表性靶基因。c 示例展示仅在 NidaWHIRLY1（NidaWHY1）、NidaMED25、NidaZIMLIKE 2（NidaZML2）中预测到的转座元件来源 CArG 盒（红色三角）。 d、e 合生心皮性状演化的潜在机制。d 黑种草、蓝花耧斗菜、中国金莲花中参与角质蜡质合成的候选基因表达谱。图示为黑种草合生雌蕊，以及蓝花耧斗菜、中国金莲花离生雌蕊。取样包含花器官第 4 至 8 发育阶段。数据为 3 个生物学重复的平均值 ± 标准差。e 3 个物种中LCR的基因结构；黑种草特有 DNA 转座子序列已高亮标注。TIR：微型反向重复转座元件（MITE）的末端反向重复序列；TSD：MITE 的靶位点重复序列。所有示例中基因结构按 5’→3’方向展示，箭头代表转录方向。源数据已作为源数据文件提供。

由于观察到的表达与序列变异并未涉及花瓣、心皮身份基因，本研究进一步聚焦其他可能参与这两项花创新性状的基因与调控通路。为解析黑种草二唇形花瓣的演化，本研究比较了黑种草、具距花瓣的蓝花耧斗菜、平瓣的中国金莲花中 AP33 同源蛋白的潜在靶基因。AP33 同源蛋白的 DNA 结合 MADS 结构域高度保守，因此本研究采用经 ChIPseq 验证的 AqcoAP33 的 CArG 盒基序，在 3 个物种基因组中扫描潜在的 AP33 结合位点。结合 CArG 盒位置信息与基因表达数据，在黑种草、蓝花耧斗菜、中国金莲花中分别鉴定到 7127、7089、4922 个候选靶基因。 黑种草的靶基因中，3826 个（53.68%）与蓝花耧斗菜（1891 个）、中国金莲花（956 个）或二者共有（979 个），这些基因主要富集于花发育、花器官数目决定、细胞分裂与生长、近轴 / 远轴极性建成等通路。剩余 3301 个为黑种草特有靶基因，显著富集在细胞周期、表皮毛乳突形成等通路。 在 NidaAP33 特有靶基因中，48.17%（1590/3301）含有转座元件起源的 CArG 盒，该比例显著高于 AqcoAP33 特有靶基因（31.46%，1085/3449）与 CochAP33 特有靶基因（37.02%，826/2231）。上述结果说明：在很大程度上依靠转座元件重塑调控网络，招募黑种草特有的 AP33 靶基因，是黑种草复杂二唇形花瓣演化的关键。

为阐明黑种草属合生心皮的演化起源（合生心皮属于后生融合，该过程可能伴随角质蜡质层破坏），本研究筛选黑种草 96 个角质蜡质合成通路候选基因，以及蓝花耧斗菜、中国金莲花中的同源基因进行比较。基因表达谱显示，黑种草中有 11 个基因在萼片、花瓣、雄蕊发育后期高表达，但在心皮中显著下调，提示其可能与心皮后生融合相关。 但其中 5 个基因在蓝花耧斗菜和 / 或中国金莲花的心皮中同样下调，说明它们并非黑种草心皮融合特有的演化因子。剩余 6 个基因中，羟基脂酰辅酶 A 脱水酶（HCD）、中链烷烃羟化酶 1（MAH1）的同源基因发生了近期串联或邻近复制；且 3 个物种里几乎所有拷贝的基因区都含有转座元件。与之不同，LACERATA（LCR）在毛茛科中维持单拷贝，且 NidaLCR 下游调控区插入了一段 DNA 转座子，而 AqcoLCR 与 CochLCR 无该特征。结合拟南芥lcr突变体可发生器官后生融合的已有证据，我们提出假说：该类群特有的转座元件插入调控区，极有可能造成 NidaLCR 在心皮特异性下调，进而推动黑种草合生心皮性状的形成。

转座元件对植物化学性状演化的贡献

β榄香烯是一种具有抗肿瘤活性的倍半萜，在毛茛科中仅富集于黑种草属的种子油中。为解析其演化机制，本研究比较了黑种草、蓝花耧斗菜、中国金莲花的倍半萜合成相关基因。结果显示，编码 β榄香烯前体合成酶的基因在 3 个物种中均存在，拷贝数变异小，且广泛表达。 与之形成对比的是，萜类合成酶 A 亚家族（TPSa）基因拷贝数差异巨大：黑种草 22 个、蓝花耧斗菜 11 个、中国金莲花仅 1 个。黑种草的两个 TPSa 基因NidaTPSa11、NidaTPSa12在富集 β榄香烯的器官（果实、种子）特异性表达；而蓝花耧斗菜、中国金莲花的 TPSa 基因在检测器官中表达量极低或检测不到。系统发育与共线性分析表明，NidaTPSa11和NidaTPSa12来自一次近期串联复制事件，该复制很可能由 LTR 转座子介导。此外，该 LTR 的插入时间估算为 3703 万年前，略早于黑种草属的分化时间。综上，LTR 介导的串联复制，叠加后续 TPSa 亚家族基因的表达分化，是黑种草 β榄香烯合成性状演化的基础。

a–d β-榄香烯演化的潜在机制。a 黑种草（Nida）、蓝花耧斗菜（Aqco）、中国金莲花（Coch）中 β-榄香烯合成相关基因的拷贝数变异。b 毛茛科 TPSa 亚家族基因的系统发育树。c 黑种草、蓝花耧斗菜、中国金莲花中 NidaTPSa11/-2/-3 周围基因组区域的共线性关系。灰色连线代表共线性同源基因。d NidaTPSa11/-2/-3 侧翼 LTR 反转座子的插入时间估算结果。粗横线代表第一、三分位数之间的四分位距（IQR），白点为中位数并标注对应数值；须线延伸至最小值与最大值。 e–h 原小檗碱类生物碱（PAs）丢失的潜在机制。e 黑种草、蓝花耧斗菜、中国金莲花中原小檗碱生物碱合成基因的拷贝数变异。f SOMT、CYP719、STOX、BBE 基因在 3 个物种中的基因组共线性。灰色连线代表共线性同源基因，目标基因对以粉色高亮；黄色、绿色连线表示 SOMT 附近一段小型倒位。g CochSOMT、CochCYP719f、CochSTOX 分别与蓝花耧斗菜、黑种草中对应的同源区域点阵图比对。黑色点代表同向序列比对，红色点代表反向序列比对。h f 图中黑种草内对应 SOMT、CYP719、STOX 区域中 LTR 的插入时间估算。粗横线代表四分位距，白点标记中位数并标注数值；须线延伸至最小值与最大值。 Mya：百万年前。源数据已作为源数据文件提供。

原小檗碱类生物碱（PAs）属于苄基异喹啉类代谢物，在黄连等毛茛科基部类群中存在，但在耧斗菜属、黑种草属等类群中缺失。为解析该现象背后的基因组基础，本研究比较了黑种草、蓝花耧斗菜与中国金莲花的原小檗碱生物碱合成通路基因。结果和预期一致：中国金莲花中全部 PAs 合成基因均存在，且在根、茎、叶、花、果实、种子中广泛表达。与之相反，另外两个物种的 PAs 合成通路并不完整：蓝花耧斗菜中关键的小檗碱桥酶基因（BBE）在所有检测器官中均不表达；黑种草中直接缺失 4 个核心基因，包括 BBE、(S)- 金黄紫堇碱 9O甲基转移酶基因（SOMT）、CYP719、(S)- 四氢原小檗碱氧化酶基因（STOX）。共线性分析显示，蓝花耧斗菜与中国金莲花在 BBE 基因座不存在共线性保守性，提示该位点发生了类群特异性染色体重排。在黑种草中，SOMT、CYP719、STOX 对应的共线性区块仍然保留，但在该区域找不到和这三个基因同源的序列，说明基因发生了完全丢失。进一步转座元件分析表明，这些基因丢失很可能发生在约 3000–3600 万年前，由 LTR 介导的异位重组造成，该时间也与黑种草属的起源时间吻合。以上结果说明，毛茛科不同类群独立丢失原小檗碱生物碱的机制并不相同，至少包括蓝花耧斗菜中 BBE 的转录沉默，以及黑种草中 LTR 介导的基因丢失。

讨论

本研究以毛茛科为研究体系，探究巨型基因组的演化过程、驱动因素及其演化后果。结果表明，毛茛科内的巨型基因组由一个约 1.46 Gb 的祖先基因组独立、逐步演化形成。黑种草基因组大小为 10.84 Gb，是目前已测序毛茛科二倍体基因组中最大的；其基因组巨型化主要源于过去 6000 万年中转座元件（TE）的持续积累，而非全基因组复制（WGD）事件。广泛发生的 TE 插入深刻重塑了黑种草的基因集：产生 4194 个携带超长内含子的基因、415 个蛋白编码序列发生显著改变的基因，以及 2473 个起源于 TE 的黑种草属特有基因。此外，TE 介导的调控改变、基因复制与丢失，很可能支撑了一系列重要特化性状的演化，包括二唇形花瓣复杂化、心皮后生融合、β榄香烯合成以及原小檗碱生物碱丢失。本研究结果为理解基因组大小演化提供了新证据，同时凸显了 TE 在驱动基因组巨型化、塑造物种特化性状中的核心作用。

毛茛科基因组大小差异可达 200 倍以上，是研究巨型基因组演化的绝佳体系。早期研究将该科染色体分为 R 型、T 型、C 型三类，并发现其和基因组大小高度相关：R 型染色体对应大基因组，T 型对应小基因组，C 型对应中等大小基因组；研究进一步推断，大型 R 型染色体与小型 T 型染色体均起源于中等大小的 C 型染色体。本研究结果不仅印证该推论，还表明毛茛科多样的基因组大小均由约 1.46 Gb 的祖先基因组演化而来。基部三个族分化后基本保留祖先基因组大小；其余支系朝着两个完全相反的方向演化：部分类群独立、逐步扩张，形成大甚至巨型基因组；而人字果族（Isopyreae）则发生剧烈收缩，形成小基因组。此外，结果显示毛茛科物种基因组大小差异主要由 TE 插入动态的差异决定，而非 TE 清除效率的不同。

毛茛科基因组大小的演化模式，和其他拥有巨型基因组的被子植物类群具有相似趋势。例如兰科基因组跨度从 341 Mb 到 54.88 Gb，推断祖先基因组大小为 4.72 Gb；后续分化中，部分支系发生基因组缩减（如拟兰亚科、石斛属、腋唇兰亚族），另一些支系逐步扩张形成巨型基因组（如兜兰亚科、鸟巢兰族、朱兰属、羊耳蒜属、布袋兰属）。与之类似，黑药花科、百合科祖先基因组均为中等大小（分别为 5.27 Gb、6.58 Gb）；演化过程中部分类群基因组大小维持稳定，另一部分发生巨型化（如黑药花科重楼族、百合科贝母属与百合属）或者基因组缩减。综上，尽管证据仍有限，但巨型基因组演化已经显现出一些共性规律：1）巨型基因组往往在不同支系中，由中小型祖先基因组独立演化形成；2）基因组大小演化通常存在分化路径，或是逐步扩张走向巨型基因组，或是剧烈收缩形成小基因组。一般认为，基因组扩张主要依靠大量 TE 积累和 / 或 WGD，为遗传创新提供原材料；而基因组缩减依靠多种 DNA 删除机制，能够提升 DNA 复制效率、降低营养消耗成本。

理解基因组巨型化的一个关键问题：它是近期 TE 爆发驱动，还是长期持续的 TE 积累驱动？越来越多证据表明 TE 是很多植物基因组的主体，TE 扩增被普遍认为是基因组巨型化的主要驱动力。已有多项研究通过完整 LTR 插入时间估算，提出很多巨型基因组来自近期 LTR 扩增，例如岷江百合（0–3 百万年前）、滇牡丹、川百合变种（0–5 百万年前）、洋葱、葱、凤丹、王百合（0–6 百万年前）。本研究采用相似策略，发现黑种草中完整 LTR 集中在近 1000 万年内积累。但祖先基因组重建结果表明，黑种草以及上述物种当前的巨型基因组，是数千万年持续扩张的结果。我们认为这种表观矛盾，主要源于仅依靠完整 LTR 进行年代估算的方法学局限，原因有两点：第一，完整 LTR 只占全部 TE 的一小部分，只能捕捉 TE 演化历史的一部分；第二，完整 LTR 容易发生降解，大多是近期产生的，仅依靠它们必然低估古老的 TE 扩增事件。

很多动物研究不局限于完整 LTR，而是对全部 TE 进行插入时间估算来解析基因组巨型化，非洲肺鱼、南极磷虾就是典型案例。这套全面的 TE 年代估算策略也已经应用于植物研究：例如云南重楼的研究发现，超过四分之三的 TE 插入发生在 2900 万年前，强烈说明古老 TE 积累对基因组扩张至关重要。将该方法应用于黑种草，我们观察到 TE 积累始于约 6000 万年前；约 92% 的 TE 插入发生在黑种草属与其姐妹类群翠雀族分化之后，并在 2700 万年前出现一次大规模扩增爆发。云南重楼与黑种草的 TE 活动动态，和祖先基因组重建得到的基因组扩张关键阶段高度吻合。此外，本研究和前人工作共同表明，不同 TE 类别演化模式不同：有的在特定时间点集中爆发快速积累，有的在漫长地质时间内平稳扩增。综上，被子植物不同支系的基因组巨型化，并非主要由单次、近期的 LTR 爆发驱动，而是数千万年间各类 TE 持续、逐步积累的结果。后续研究探索这种持续 TE 增殖的诱因（古气候事件、历史种群动态或其他未知因素），将进一步加深我们对巨型基因组演化的理解。

黑种草基因组组装结果，加深了我们对 TE 塑造基因与性状演化的认知。虽然只有一小部分 TE 位于基因区域，但它们深刻影响蛋白编码基因的演化。首先，TE 显著增大内含子长度，黑种草内含子平均长度达到 3016 bp，并产生 4194 个富集 TE 的超长内含子基因。之前百合、油松的研究显示，超长内含子不一定会损害转录效率。在黑种草中我们也观察到类似现象，并且发现 46.38% 携带富集 TE 超长内含子的基因表达模式保守，该比例高于内含子 TE 较短的基因。这种高转录活性与表达保守性，可能和富含 TE 的内含子中更高的 DNA 甲基化水平有关，油松、滇牡丹、王百合中也报道过类似现象。其次，TE 改变蛋白编码序列、重塑基因组成。具体来说，黑种草属中发生的 TE 插入，通过起始 / 终止密码子改变、外显子插入、外显子化、内含子保留等多种机制，使 415 个原有基因编码蛋白多出一段氨基酸序列；TE 插入还产生了 2473 个大部分序列起源于 TE 的黑种草属特有基因。这类谱系特有基因的功能仍有待深入解析，但超过半数能够表达，部分存在组织特异性，说明它们在黑种草中具备潜在功能。

TE 介导的顺式调控分化推动花瓣复杂化与心皮融合。黑种草花瓣为二唇形，带有长绒毛、短表皮毛、假蜜腺等高度特化附属结构，和蓝花耧斗菜的距状花瓣、中国金莲花的扁平花瓣明显区分；但这三类花瓣都起源于雄蕊，也被称为蜜腺叶。已有研究认为毛茛科花瓣形态差异主要来自近轴 / 远轴基因表达变化、区域生长分化以及基因共选择。本研究发现，虽然大约一半 AP33 靶基因在物种间共享，但物种特有靶基因和各花瓣类型独有的特化性状相关。例如黑种草 NidaAP33 特有靶基因包含参与二唇结构形成（如NidaAS2）、表皮毛发育（如NidaCPC）、生长调控（如NidaARF8、NidaMED25）、色素合成（如NidaLUT2）的基因；蓝花耧斗菜 AqcoAP33 特有靶基因包含调控距发育（如AqcoPOP）、叶绿素降解（如AqcoPPH、AqcoPAO、AqcoACD2）的基因。更重要的是，48.17% 的 NidaAP33 特有靶基因带有 TE 起源的 CArG 盒，该比例显著高于另外两个物种。此外，本研究结果显示，NidaLCR在心皮特异性下调（很可能由 DNA 转座子介导调控区改变引起），可能促成合生心皮性状的起源。上述结果展示了 TE 如何重连基因调控网络、驱动黑种草形态性状演化，当然还需要后续功能验证。

TE 在植物化学性状演化中扮演双刃剑角色：一方面通过基因复制催生 β榄香烯这一新性状，另一方面通过基因缺失造成原小檗碱生物碱丢失。已有观点认为，和次生代谢相关的串联重复基因，通常产生于不等交换；不等交换在 TE 富集区域、TE 含量超过 70% 的基因组中高频发生。本研究结果支持该观点：黑种草中 β榄香烯合成候选基因 NidaTPSa11、NidaTPSa12 为串联复制，且位于 LTR 高度富集的区域。TE 除了助力新植物化学性状获得，TE 介导的不等交换同样会造成性状丢失。具体而言，黑种草中 SOMT、CYP719、STOX 等原小檗碱生物碱合成基因的丢失，很可能由 LTR 介导的异位重组导致。综上，TE 在驱动基因与性状演化方面具备巨大潜力，后续仍需要更多研究解析 TE 在巨型基因组中更加广泛的功能图谱。

方法

植物材料与培养条件

黑种草（N. damascena）与蓝花耧斗菜（Aq. coerulea）种子播种于营养土中，培养条件：温度 25 ℃、相对湿度 60%、光周期 12 h 光照 / 12 h 黑暗。蓝花耧斗菜幼苗长至 4–6 片真叶时，于 4 ℃短日照（8 h 光照 / 16 h 黑暗）条件下春化一个月，之后放回初始培养条件。中国金莲花（C. chinensis）购自湖北利川。

祖先基因组大小重建

为解析毛茛科巨型基因组的演化，本研究首先利用流式细胞术测定 20 种毛茛科植物的基因组大小（附图、补充数据）。对每个物种，搭配已知基因组大小的合适内参，取约 30 mg 新鲜叶片，置于 1 mL 冰浴 Galbraith 缓冲液（45 mM MgCl₂、20 mM MOPS、30 mM 柠檬酸钠、0.1%（v/v）Triton X100，pH 7.0）中切碎；所得细胞悬液经 400 目尼龙滤膜过滤，加入碘化丙啶（终浓度 100 μg/mL）染色 10 min。染色后的细胞核采用 MoFloXDP 流式细胞仪（美国 Beckman Coulter）检测，使用 Summit v5.2 软件计算基因组大小。另外，从已发表文献中收集 30 种毛茛科物种以及 18 个外类群物种的基因组大小数据（补充数据，内含参考文献）。

构建包含上述物种的时间标定系统发育树（附图），构建思路参照 Timetree 项目。简要来说：科、族水平的拓扑结构与分化时间主干树采用 He 等人的综合分析结果；再利用 Mesquite v3.70 软件的树编辑工具，将已发表、可信度高的族内 / 属内系统发育分支手动嫁接至主干树上。随后在 Mesquite 中调整分支长度，匹配原始文献报道的分化时间，完成时间标定。外类群选自毛茛目下 6 个非毛茛科的科，均具备可靠基因组大小数据、解析清晰的系统发育关系。

为筛选适合祖先基因组大小重建的演化模型，使用 R 包geiger v2.0.11 的fitContinuous函数，评估 4 种常用模型：布朗运动模型（BM）、奥恩斯坦乌伦贝克模型（OU）、早期爆发模型（EB）、白噪声模型（WN，无系统发育信号，作为对照）。依据加权校正赤池信息准则（AICc），选择 BM 模型为最优模型（补充数据）。基于该模型，使用 R 包phytools v1.9.16 的anc.ml函数推断各演化节点的祖先基因组大小，并通过contMap可视化。

基因组测序、组装与质量评估

黑种草基因组测序：取幼嫩叶片，使用天根 DNAsecure 植物基因组 DNA 提取试剂盒提取基因组 DNA。采用 Illumina NovaSeq 平台获得 1018.20 Gb 的 150 bp 双端短读长数据；PacBio Sequel II 平台获得 274.82 Gb 单分子长读长数据（补充数据）。Hi-C 测序：幼苗经甲醛固定，交联 DNA 使用 MboI 限制性内切酶酶切；Illumina NovaSeq 150 bp 双端测序，共得到 247.04 Gb 有效 Hi-C 数据（补充数据）。

分别采用流式细胞术（如上所述）与 kmer 分析估算黑种草基因组大小。kmer 分析使用 Jellyfish v2.2.10 统计 Illumina 短读长的 21mer 频率分布，再用 GenomeScope v2.0 预估基因组大小。基因组组装：PacBio HiFi 读长采用 hifiasm v0.2.0 默认参数从头组装；使用 Purge_dups v1.2.5 去除冗余序列。基于 Hi-C 互作图谱，使用 LACHESIS（参数CLUSTER_N = 6）将 contig 锚定并定向到 6 条假染色体上；借助 Juicebox v1.9.0，依据 Hi-C 互作频率手动校正组装错误的 contig。

组装质量评估包含三项：1）BWA v0.7.17（参数mem -M -k 30）将 Illumina 双端读长回贴到组装基因组，评估一致性；2）BUSCO v5.4.7，使用embryophyta_odb10数据库评估基因组完整度；3）LAI 评估组装连续性。

基因组注释

采用三种策略预测黑种草蛋白编码基因：

三种来源的基因模型通过 EVidenceModeler（EVM）整合。整合后的基因集通过 BLAST（E 值≤1×10⁻⁵）比对多个数据库完成功能注释：SwissProt、NT、NR、Pfam、eggNOG、GO、KEGG。转录因子利用 PlantTFDB v5.0 鉴定。

重复序列注释分为两套策略：

两套重复库合并，使用 RepeatMasker v4.1.2-p1 注释转座元件（TE）。LTR_retriever v3.0.4 识别完整 LTR 与 solo-LTR；solo-LTR 再经 TEsorter v1.5.1 进一步过滤；利用 BEDTools v2.31.0 从全部 LTR 中剔除完整 LTR 与 solo-LTR，得到截断型 LTR。tRNA 使用 tRNAscan-SE v1.3.1 默认参数预测；rRNA、snRNA、miRNA 通过 Infernal v1.1.3（默认参数）比对 Rfam 数据库鉴定。

全基因组复制（WGD）事件鉴定

为鉴定黑种草的 WGD 事件，开展基因组内共线性分析，同时与另外 6 种毛茛目物种基因组进行跨物种共线性比较。蛋白序列使用 BLASTP v2.5.0+（E 值≤1×10⁻⁵）进行全对全比对，得到候选同源基因；BLAST 结果输入 MCScanX（默认参数）识别共线性区块。对每个共线性区块，使用 KaKs_Calculator v2.0 的 NeiGojobori 算法计算锚定基因的 Ks；取锚定基因的中位数 Ks 表征该区块分化水平。

WGD 发生时间估算的标尺：黑种草与虞美人分化时间为 109.6 百万年，对应的同源 Ks 峰值为 0.92。参照 Wang 等人的方法对 Ks 峰值进行校正： （1）利用蓝花耧斗菜旁系同源 Ks 峰值，统一校正所有物种的毛茛目古老 WGD 事件； （2）利用虞美人与R. cassubicifolius的直系同源 Ks 峰值，校正罂粟科与其他毛茛目支系的分化； （3）利用蓝花耧斗菜和黑种草的直系同源 Ks 峰值，校正蓝花耧斗菜–R. cassubicifolius、黑种草–R. cassubicifolius的分化时间。

转座子插入时间估算

采用两种策略估算转座子（TE）的插入时间。第一种策略专门针对完整 LTR 转座子：比对其 5' 端和 3' 端末端重复序列，基于 Kimura 双参数（K2P）模型计算核苷酸距离（K）。第二种策略适用于全部已注释 TE：将每条 TE 序列与其对应的一致性序列（consensus sequence）比对，计算 K 值。两种策略下，TE 插入时间（T）均按下式计算： （1） 式中，替换速率 r = 2.85 × 10⁻⁹，即每个位点每年的平均替换率。r 值计算：选取 5 种毛茛科物种、2 个外类群（虞美人、E. pubescens）共 1007 个单拷贝基因的核苷酸序列，使用 r8s v1.81 默认参数运算；以罂粟科–毛茛科分化时间（102.9–117.2 百万年前）、小檗科–毛茛科分化时间（84.0–97.1 百万年前）作为校正标定点。

直系同源基因鉴定与基因获得 / 丢失分析

从已发表文献获取 33 个物种的蛋白编码序列（30 种毛茛科、2 种小檗科、1 种罂粟科，补充数据含参考文献）。将上述序列与黑种草序列合并，使用 OrthoFinder v2.5.4（参数-S blast -M msa -T raxml -t 128）聚类得到基因家族 / 亚家族（同源组 orthogroup）。基于同源组聚类结果，使用 Count 软件的 Wagner 简约法推断基因家族获得与丢失事件。

基因表达谱分析

收集黑种草、蓝花耧斗菜、中国金莲花成熟器官（根、茎、叶、花、果实、种子）以及不同发育阶段花器官的 RNA-seq 数据，每个物种、每个发育阶段、每个器官设置 3 个生物学重复。黑种草与中国金莲花的花发育阶段参照蓝花耧斗菜对应的形态特征划分。黑种草花器官、中国金莲花根 / 茎 / 叶的 RNA-seq 数据取自已发表文献；其余 174 个样本为本研究测定（补充数据）。每个样本按前述方法提取总 RNA 并测序，平均得到 6 Gb 的 150 bp 双端干净 reads。干净 reads 使用 HISAT2 v2.21 默认参数比对到对应参考基因组；采用 Stringtie v2.21 计算每个样本基因的 TPM（每百万转录本），取 3 个重复的平均值代表基因表达量。

为实现跨物种、跨器官的表达量可比，参照 Brawand 等人的方法对平均 TPM 做标准化：首先筛选黑种草、蓝花耧斗菜、中国金莲花之间严格一对一的高度保守单拷贝基因；依据各样本中这些基因的 TPM 排序，筛选 147 个排序最稳定的基因（在所有样本中排名稳定落在四分位距内，补充数据）；计算每个样本中该组基因的中位数 TPM，推导缩放系数，统一校正到基准值；最后用该缩放系数校正样本内全部基因的 TPM。除非另行说明，本文展示的所有 TPM 均为校正后数值。

转座子插入对基因的影响评估

依据 TE 与基因组注释区域（CDS、内含子、侧翼序列）的重叠关系定义 TE 插入位点，统计规则：若一条 TE 完全落在单个区域内，记为 1 次插入，TE 全长作为插入长度；若一条 TE 跨多个区域（如同时覆盖 CDS 和侧翼区），在每个重叠区域各记 1 次插入，记录落在该区域内的片段长度。

选取黑种草、蓝花耧斗菜、中国金莲花 10 组匹配器官（6 种成熟器官：根、茎、叶、花、果实、种子；4 类花器官：萼片、花瓣、雄蕊、心皮）的表达谱，评估内含子 TE 积累对基因表达的影响。萼片、花瓣、雄蕊、心皮取该器官全部发育阶段里最高 TPM 作为表达水平。基因表达保守的判定标准：黑种草该基因在所有器官的表达模式和另外两个物种的直系同源基因相近；若某基因仅在某一个器官相比同源基因表达差异≥2 倍、其余器官表达保守，则判定该基因在黑种草上调或下调；器官特异性表达改变：相比同源基因，该基因在部分器官显著上调、另一些器官显著下调。使用 R v4.3.2 进行卡方检验，分析内含子 TE 含量与基因表达分化的关联。

两步鉴定黑种草属特有 TE 来源基因：①筛选黑种草中蛋白编码序列部分或全部与 TE 重叠的基因；②按前述 OrthoFinder 流程查找其直系同源基因，仅在黑种草属以外找不到同源基因的，定义为黑种草属特有 TE 衍生基因。

花性状与植物化学性状相关基因的系统发育分析

以黑种草花发育 MADS-box 基因、拟南芥角质蜡质与倍半萜合成基因、中国金莲花原花青素（PA）合成基因为查询序列，BLASTP v2.5.0+（E≤1×10⁻⁵）检索毛茛目物种基因组 / 转录组中的同源序列。每个基因家族的氨基酸序列用 MAFFT v7.505 比对，再用 PAL2NAL v14 生成基于密码子的核苷酸比对矩阵；IQ-TREE 2.1.4-beta 构建最大似然树（参数-bb 1000 -nt AUTO）。

AP3-3 直系同源基因潜在靶基因预测

使用 AlphaFold Server 预测 NidaAP3-3、AqcoAP3-3、CochAP3-3 蛋白结构。采用 4 个经 ChIP-seq 验证的 AqcoAP3-3 结合 CArG 基序，扫描黑种草、蓝花耧斗菜、中国金莲花基因组，寻找 NidaAP3-3 与 CochAP3-3 的候选靶基因。候选靶基因定义：调控区（翻译起始位点上游 3 kb 至翻译终止位点下游 1 kb，排除 CDS）至少包含 1 个 CArG 盒，且在花瓣中表达。若 CArG 盒与 TE 序列部分或完全重叠，则该 CArG 盒判定为 TE 起源（附图）。

GO 富集分析

对目标基因集做 GO 富集时，使用其拟南芥同源基因作为输入；采用 g:Profiler 进行富集，Benjamini–Hochberg 校正 FDR 阈值 0.05，保留显著富集的生物学过程条目。

Data availability

The genome sequencing data for N. damascena and the transcriptome sequencing data for N. damascena, Aq. coerulea, and C. chinensis have been deposited to National Genomics Data Center under accession PRJCA041492. The genome assembly and annotations of N. damascena are available at Figshare [https://figshare.com/s/7e80a04b8c87b97b24ea]. Source data are provided with this paper.