红花变豆菜叶绿体基因组组装与序列特征分析研究

红花变豆菜Sanicula rubriflora Fr. Schmidt是伞形科（Apiaceae）变豆菜属多年生草本植物，是一种可食用的山野菜，具有一定的药用价值，变豆菜属全世界约40余种[1]。其生长在山间林下，阴湿及腐殖质较多的地方，海拔200～470 m，主要产地在我国的东北地区，另外，蒙古、朝鲜、日本北部也有分布[2]。变豆菜属植物在中国约有18种，具有一定的经济价值和药用价值[3]，薄片变豆菜S. lamelligera Hance又叫血经草，具有治疗跌打损伤和风寒感冒等功效[4]，天蓝变豆菜S. caerulescensFranch.、直刺变豆菜S. orthacantha S. Moor和川滇变豆菜S. astrantiifoliaWolff ex Kretsch均为有名的地方药材，具有活血化瘀和止咳化痰的药效，同属的红花变豆菜也可能会具有相应的药用价值[5]。变豆菜属植物主要基于传统的形态学进行鉴定，如从花粉形态[6]和果实微形态出发进行分类，但是部分品种只从形态学鉴定很难区分，导致其传统分类发展缓慢，因此需要寻找新的方法来对其鉴定手段进行补充。

叶绿体普遍存在于陆地植物、藻类和部分原生生物当中，是绿色植物进行光合作用的细胞器[7]，功能包括产生色素、合成糖和某些氨基酸等[8]。被子植物中，叶绿体基因组相对保守，绝大多数为双链环状结构，包括1个小的单拷贝区（small single copy，SSC）、1个大的单拷贝区（large single copy，LSC）以及2个编码相同，方向相反的序列（inverted repeats sequence，IR）IRa和IRb [9]。相对于核基因组，叶绿体基因组具有结构保守、碱基变异速率适中、易于测序等多种优势，已广泛应用于各植物类群的系统进化研究[10]。基于本课题组前期对红花变豆菜叶绿体基因组的研究[11]，对其序列数据进一步挖掘分析，对组装后的叶绿体基因组结构、长重复序列分析、简单重复序列（simple sequence repeat，SSR）位点、密码子偏好性以及与其他伞形科植物使用最大似然法（maximum likelihood，ML）构建系统发育树，分析近缘物种亲缘关系，为红花变豆菜的鉴定、系统发育研究提供新的方法。

1材料与方法

1.1 材料

新鲜的红花变豆菜样本于2020年8月20日采集自中国黑龙江省伊春市（N47°81′08′′，E128°90′97′′）。样品保存在黑龙江中医药大学，样品标本号为YCL20190507007，由黑龙江中医药大学马伟研究员鉴定为伞形科变豆菜属植物红花变豆菜S. rubrifloraFr. Schmidt。将采集的红花变豆菜新鲜幼嫩的叶片，通过液氮速冻，并进行研磨，于−80 ℃冰箱保存备用。

1.2方法

1.2.1DNA 提取使用CTAB法提取红花变豆菜的全基因组DNA。

1.2.2测序总DNA样品经武汉贝纳科技服务有限公司检测合格后，用机械打断的方法（超声波）将DNA片段化，然后对片段化的DNA进行片段纯化、末端修复、3′端加A、连接测序接头，使用琼脂糖凝胶电泳进行片段大小选择，进行PCR扩增形成测序文（NEBNext ® Ultra™ DNA Library Prep Kit for Illumina ® ），建好的文库先进行文库质检，质检合格的文库用Illumina NovaSeq进行测序，对测序得到的序列，通过SOAPnuke（version：1.3.0）软件进行低质量数据过滤，过滤标准为去除N碱基含量超过5%的reads、去除低质量（质量值≤5）碱基数目达到50%的reads和去除有adapter污染的reads，最终得到有效数据（clean reads）。

1.2.3叶绿体基因组组装使用SPAdes（version：3.13.0；参数：−k 127）[12]软件进行基因组拼接，选用NCBI数据库中近缘参考植物变豆菜叶绿体基因组（MK208987.1）与拼接结果进blastn（version：BLAST 2.2.30＋；参数：−evalue 1×10 −5 ）比对，后续使用Gapcloser（Version：1.12）对gap（含N序列）进行补洞。

1.2.4叶绿体基因组特征分析利用专门针对叶绿体的注释软件CPGAVAS2（http://47.96.249.172: 16019/analyzer/annotate）[13]进行基因注释并绘图红花变豆菜叶绿体基因组物理图；使用VMATCH（http://www.vmatch.de/）软件（参数：minimal repeat size 30 bp）查找叶绿体基因组中的散在长重复序列片段；采用MISA软件[版本：1.0；默认参数；对应的各个重复单元（unit size）的最少重复次数分别为1-8、2-4、3-4、4-3、5-3、6-3]对叶绿体进行SSR检测；使用CodonW（Version：1.4.4）对密码子偏好性进行分析，统计估算相对同义密码子的使用频率；通过NCBI数据库下载6种伞形科植物，包括变豆菜Sanicula chinensisL.（MK208987）、明日叶Angelica keiskeiL.（MW125613）、隔山香Ostericum citriodorumL.（MT501096）、辽藁本Ligusticum jeholenseL.（MN652885）、防风Saposhnikovia divaricataL.（MN857472）、白芷Angelica dahuricaL.（KT963037）以及2种非伞形科植物玉米Zea maysL.（NC001666）和拟南芥Arabidopsis thalianaL.（NC000932）的叶绿体基因组序列，使用MAFFT [14]软件进行多序列比对，最后使用MEGA X [14]软件基于ML法构建系统发育树。

2结果与分析

2.1红花变豆菜叶绿体基因组测序及结构解析

基于课题组已报道的研究发现，红花变豆菜的叶绿体基因组大小为155 700 bp，具有典型的环状四分体结构：包括1个LCS区（85 979 bp），1个SCC区（17 053 bp）和2个反向互补的IRs区（26 333 bp），GC含量为61.83% [14]。

后续通过对红花变豆菜叶绿体基因组物理图谱的绘制（图1）及进一步统计分析发现，测序有效数据Q20值为97.41%、Q30值为92.84%，根据比对覆盖度发现变豆菜（MK208987.1）为最优参考序列，其覆盖度为96.44%，不含有gap。IRs、LSC和SSC 4个区域的GC值都存在一定的差异，IRs区域的CG值最高，为42.92%；SSC区域的CG值最低，为32.56%，LCS区域的CG值介于二者之间，为36.39%（表1）。

2.2 叶绿体基因组结构基本特征分析

基于课题组已报道的研究发现，红花变豆菜叶绿体基因组共成功注释了130个基因，其中蛋白编码的基因86个，tRNA基因36个，rRNA基因8个[11]。

通过对以上结果进一步的统计分析发现，根据其基因编码功能可以分为4类（表2），第1类是与自我复制有关的基因，分为5个组，分别是tRNA基因、核糖体小亚基、核糖体大亚基、RNA聚合酶亚基和核糖体RNA基因，成员数量分别为26、12、9、4和4个；第2类为与光合作用相关的基因，细分为7组，分别是光系统I、光系统II、ATP合成亚基、细胞色素亚基、二磷酸核酮糖羧化酶大亚基、NADPH乳酸脱氢酶亚基和依赖ATP蛋白酶P基因亚基，成员数量分别为5、15、6、6、1、11和1个；第3类为其他功能基因，包括成熟酶、外膜蛋白基因、C型细胞色素合成基因和乙酰辅酶A羧化酶亚基各1个；第4类为4个未知功能的基因，还需要进一步研究，以确定其功能。

2.3 长重复序列分析

所有的长重复序列包括3种类型，包括正向重复序列（forward repeat sequence）、回文重复序列（palindromic repeat sequence）和串联重复序列（tandem repeat sequence），这些长重复序列可能具有促进叶绿体基因组重排的功能，并且可以增加其居群遗传多样性[16]。在红花变豆菜叶绿体基因组中共发现了1个正向重复序列，3个回文重复序列，并未发现串联重复序列。其中1个正向重复序列和2个回文重复序列非常短，分别为34、34、30 bp，1个回文重复序列非常长，为26 333 bp，是叶绿体基因组的IR区域。利用这些重复序列，可以为以后开发种群进化标记研究提供基础。

2.4 SSR位点分析

SSR又称短串联重复标记、微卫星序列标记[17-18]。SSR表示由基因组中1～6个核苷酸组成的基本单位重复多次构成的一段DNA序列，广泛存在于基因组的各个区域，且侧翼序列通常都是保守性较强的单一序列，可用于个体或物种之间的多态性研究[19]。通过对红花变豆菜叶绿体基因组的分析，共发现了168个SSR位点，分为4个类型，包括复杂重复类型112个，其数量最多；3个碱基重复类型2个，其数量最少；2个碱基重复类型47个；4个碱基重复类型7个，没有单碱基重复类型。97.60%的SSR位点都含有A/T的碱基，仅有4个SSR位点由G/C碱基组成，说明A/T碱基具有碱基偏好性，这可能跟红花变豆菜叶绿体基因组中A/T碱基含量占比高（61.82%）有关联，造成这种偏好性的原因可能与解链难易程度有关。

2.5 密码子偏好性分析

红花变豆菜叶绿体基因组中86个蛋白由21 963个密码子共同编码（表3）。由AUU编码的异亮氨酸数量最多，为876个；由UGC编码的半光氨酸数量最少，仅有60个。3种终止密码子UAA、UAG和UGA在红花变豆菜的叶绿体基因组中使用，其在同义密码子中的占比为50%、25%和25%。通过对86个蛋白编码基因序列分析，得出了同义密码子相对使用度（relative synonymous codon usage，RSCU）[20]，RSCU是对同义密码子使用偏好的评估，通过统计发现，由UUA编码的亮氨酸使用最频繁，偏好性最大；由AGC编码的丝氨酸使用频率最低，偏好性最小；由UGG编码的色氨酸和由AUG编码的甲硫氨酸无偏好性。

2.6 红花变豆菜与其他7种伞形科植物叶绿体基因组的比较分析

本研究基于近些年报道的直刺变豆菜[21]、变豆菜[22]、明日叶[23]、隔山香[24]、辽藁本[25]、防风[26]和白芷[27]共7种伞形科植物与红花变豆菜的叶绿体序列长度、编码基因数量及组成结构进行比较（表4）。经过数据统计分析，伞形科内8种植物叶绿体基因组的长度为146 918～155 919 bp，并且都是由典型的四分体结构组成；其大拷贝区和小拷贝区的序列长度差异不大；除了辽藁本的反向重复序列长度为明日叶和防风的2倍，其他7种反向重复序列长度差异也不大；从编码的蛋白质数量来看，红花变豆菜编码的蛋白质数量最多，可能其叶绿体在生长过程中能行使更多的功能。红花变豆菜相比与变豆菜与直刺变豆菜相比，编码的蛋白最多和SSC区序列最长。除此之外，这8种伞形科植物的叶绿体基因组基因顺序和结构与大多数已报道的被子植物叶绿体基因组相似，这说明在植物进化过程中，叶绿体基因序列具有高度的保守性[28]。

2.7叶绿体基因组全序列的聚类分析

叶绿体基因组的聚类分析对植物发育进化研究具有重要意义[28]。通过对红花变豆菜与其他7种伞形科植物叶绿体基因组序列长度，基因结构以及编码蛋白分析，发现其并未有太大差异，于是通过ML法对其构建系统发育树，并加入玉米和拟南芥的叶绿体基因组序列作为外类群，对红花变豆菜进行进化分析，以确定其在伞形科植物中的进化位置（图2）。本课题组发现不同科的植物分为不同的进化分支，8种伞形科植物聚在了1个分支上；说明植物科间关系明确；在这个大的分支基础上，红花变豆菜与变豆菜和直刺变豆菜聚到了1个小分支上，并且自展值为100，说明这3个物种在漫长的进化过程中亲缘关系最为接近。

3 讨论

本研究采用高通量测序技术对红花变豆菜叶绿体基因组进行重测序，并以已报道植物变豆菜叶绿体基因组为参考，成功组装出其完整的叶绿体基因组。相对于传统意义上叶绿体基因组序列的获取，如紫荆泽兰Ageratina adenophoraL. [29]和菝葜Smilax china L. [30]都是采用先从植物样本中分离出叶绿体，然后对其叶绿体进行DNA提取，最终再通过测序技术，实现叶绿体基因组序列的获取，该方法复杂繁琐、不易操作，不利于大范围使用。利用高通量测序技术对红花变豆菜全基因组DNA进行重测序，省略了先分离植物叶绿体再提取叶绿体DNA的复杂操作过程[31]，只需提取红花变豆菜全基因组DNA，进行高通量测序，选取已报道植物变豆菜叶绿体基因组序列作为参考基因组，将其所测得的全基因组序列与参考叶绿体基因组序列进行BLASTN比对，提取出关联的叶绿体raw reads，使用过滤软件SOAPnuke对reads进行低质量序列过滤得到clean reads，再利用SPAdes软件对这些序列进行组装及优化，最后使用GapCloser软件对组装结果进行补洞，最终得到完整的红花变豆菜叶绿体基因组序列。本研究使用的方法相较于传统方法，大大简化了试验步骤、降低了试验所需的时间、减少了实验成本，并且该方法实用性广，限制条件少，为大量植物测序叶绿体基因组提供了可能。

变豆菜属Sanicula L.在克朗奎斯特系统（Cronquist system）、哈钦松系统、恩格勒系统（Engler system）、被子植物发育系统（angiosperm phylogeny group，APG）以及《中国植物志》等各个国家的植物分类体系中都被分为一个自然的类群，归属于伞形科（Apiaceae）、变豆菜亚科（Saniculoideae Drude）[32]。其中一些变豆菜属植物具有一定的药用价值，如变豆菜具有治疗风湿、咳嗽和激活血液循环的功效，直刺变豆菜具有祛风止咳、活血通络和清热解毒的功效，也是一味著名的民族药，收录于《四川省中草药标准(试行稿)》。中国变豆菜属植物约有18种，属于特有的为12种，但是鉴定困难，种间属间易混淆，极其不利于此类植物的发展与利用。针对其分类已经有相关的研究，包括通过形态学研究和分子系统学研究[34]。本研究通过对红花变豆菜叶绿体基本组进行测序及数据挖掘，分析了红花变豆菜的长重复序列、SSR位点和密码子偏好性，为变豆菜属植物研究提供了更为详细、完善的资料, 为该物种叶绿体基因工程和系统进化分析提供参考依据。