日本落叶松SNP和InDel位点及其所在抗生物胁迫基因挖掘

摘要：[目的]挖掘日本落叶松SNP和InDel位点及其所在抗生物胁迫基因，为日本落叶松分子育种提供分子标记和候选基因。[方法]使用158份来自中国、日本和英国地区的日本落叶松转录组数据，首先进行分子标记位点的鉴定和分类，随后比较了活动期和休眠期基因的表达水平，最后将带有可靠分子标记的差异表达基因进行功能注释。[结果]本研究共鉴定到515 935个SNP位点和1 056个InDel位点，它们分布在35 827个基因上。通过比较不同地区的位点数量，推测日本地区的日本落叶松遗传多样性较为丰富。至少在50份转录组中出现的非同义突变SNP位点有6 444个，InDel位点为10个，它们分布在3 742个基因上，可以作为可靠的分子标记进行利用。活动期和休眠期的转录组比较后，发现带有可靠分子标记的2 569个基因差异表达；GO注释后发现其中101个基因与植物对真菌、细菌、卵菌、病毒、昆虫和线虫的抗性反应有关。[结论]这些结果不仅为日本落叶松全基因组关联分析和全基因组选择育种提供了分子标记，也为利用转基因和基因编辑手段进行遗传育种提供了候选基因。

关键词：日本落叶松；分子标记；抗性基因；转录组

中图分类号：S722;S791.223 文献标识码：A 文章编号：1001-1498（2024）05-0054-11

单核苷酸多态性标记（single nucleotidepolymorphism，SNP）是由单碱基的转换及颠倒等突变而产生的不同个体之间碱基变异，是物种可遗传变异中最常见的一种，目前已被应用于植物抗性育种研究中，是全基因组关联分析的基础。Calic等人在山毛榉（Fagus grandifolia Ehrh.）中检测出4个与抗病性状极显著相关的SNP。插入缺失性标记（insertion and deletion，InDel）是指同一物种不同个体之间基因组同一位点的序列发生不同大小核苷酸片段的插入或缺失，目前在植物抗性育种研究中报道较少。

日本落叶松（Larix kaempferi （Lamb.） Carr.）原产日本，被引种到世界各地作为重要造林树种。目前，关于日本落叶松林染病的案例越来越多，例如因松材线虫（Bursaphelenchus xylophilus）引起的松材线虫病，因疫霉（Phytophthoraramorum）引发的疫霉病。为了适应生存环境，在一个生长周期中日本落叶松需要经历两个时期：活动期和休眠期。在活动期，落叶松生长旺盛，能够为病原物提供丰富的资源和营养，且生长环境也更适宜病原物生长和繁殖；而在休眠期，落4p/YAJ1nS+ij3rSQpcwiEw==叶松停止生长，为病原物提供的资源和营养有限，再加上环境的影响，病原物的活动和繁殖受到限制。在这两个时期中有很多基因差异表达，构成了日本落叶松适应生存环境的分子基础，因此“在活动期和休眠期差异表达的部分基因是否在日本落叶松防御病原物攻击的过程中发挥重要作用”以及“这些基因中存在的序列变异是否可能造成不同基因型拥有不同抗性”值得探讨。

日本落叶松基因组序列[9-10]以及大量转录组的测定，为鉴定序列变异提供了丰富的序列资源，也为推动分子标记辅助育种技术在日本落叶松抗性育种的应用奠定基础。本研究利用日本落叶松158份转录组数据对SNP和InDel位点进行鉴定和分析，并比较了这些位点所在基因在活动期和休眠期的表达，以期为日本落叶松抗性育种提供分子标记和候选基因。

1 材料和方法

本研究的全流程如图1所示。

1.1 转录组数据收集

78份转录组测序材料采自中国林业科学研究院，液氮冷冻后-80℃超低温冰箱保存，干冰运至北京诺禾致源科技股份有限公司进行转录组测序，转录组数据已提交到国家生物信息中心（https：//www.cncb.ac.cn）数据库（登录号：PRJCA022796）。并且在公共数据库NCBI（www.ncbi.nlm.nih.gov）上，筛选和下载测序平台尽量一致、Strategy为RNA-Seq、Layout为PAIRED的日本落叶松转录组数据（表1），共同用于鉴定分子标记。参考转录本来自日本落叶松的50 690条编码序列（Codingsequence， CDS）（http：//btg.kazusa.or.jp/blast.html）。

1.2 转录组数据处理

利用fastp软件（v0.32）进行数据过滤，去掉测序接头，保留长度在20 bp以上的序列；将过滤后的数据按照长度进行筛选，去掉长度小于30bp，或者只有一端的read。将质控后得到的cleanreads，使用BWA软件比对到参考转录本上，获取reads在参考转录本上的定位信息。使用Salmon软件统计比对到每个转录本上的read数量，并计算每个转录本的表达量，使用Transcripts perkilobase million（TPM）表示。

1.3 基因组变异检测和分类

使用变异检测软件bcftools，从1 58份转录组数据中分别检测SNP和InDel。随后对原始结果进行过滤。其中mpileup设定参数：-C50-d 80 000 -L 10 000-q 20；call设定参数-vMm；filter设定参数-i＇%QUAL＞10＆（DP4[2]+DP4[3]）＞2‘。最后使用SnpE仟工具对检测到的变异进行功能分类。使用TBtools软件，对来自不同地区（中国黑龙江，中国辽宁，日本和英国）转录组中的SNP和InDel进行韦恩分析。

1.4 差异基因分析

158份转录组数据包含本课题组之前测定的12份来自活动期和休眠期的日本落叶松转录组数据。本研究使用TBtools软件对这12份转录组数据进行了差异基因分析（p＜0.05），随后利用eggNOG和SwissProt数据库对得到的差异基因进行注释。

2 结果

2.1 转录组数据统计

158份转录组数据统计显示，每份数据大小在17 292 864～159 873 888 bp之间；质控后得到的数据大小在16 053 498～155 706 686 bp之间，Q30在80.76%～97.3g%之间，高质量数据比例在80.89%～99.74%之间。将质控后得到clean reads Lk对到参考转录本上，发现共有48 094个转录本得到比对，占总转录本数量的94.88%，每份转录组中有55.72%～88.00%的clean readsLC对到参考转录本上。

2.2 SNP、InDel鉴定及分类

每份转录组中SNP的数量在50 520～154 343个之间，InDel的数量在3 045～9 457个之间。158份转录组共有SNP位点515 935个（表2），共有InDel位点1 056个（表3），这些位点分布在35 827个转录本上，占总转录本数量的70.68%。在SNP中，注释为错义突变的有302 087个，占比最大，为58.55%，注释为同义突变的有203 529个，占比39.45%；注释为终止密码子获得的有9 091个，注释为起始密码子丢失的有318个；同时注释为两种类型的有861个。在InDel中，注释为移码突变的最多，有742个，占比70.27%；同时注释为两种类型的有34个，且其中一种类型都为终止密码子获得。

分别从日本、英国、中国辽宁和中国黑龙江来源的转录组中随机选择6份进行SNP和InDel位点分析。其中，日本独有的位点有52 798个，英国独有的有25 785个，中国黑龙江独有的有46 1 70个，中国辽宁独有的有24 056个，这4个地区共有的有12 879个（图2）。

2.3 位点筛选

CDS内的非同义突变可能会导致基因功能的改变。本研究筛选出注释为非同义突变的SNP位点有312 399个，InDel位点全部为非同义突变，有1 056个，这些位点所对应的基因共有33 101个。在158份转录组中都出现的位点极少，SNP位点仅有755个（0.24%），InDel位点仅有1个（0.09%）。在158份转录组中，至少在25份里出现的SNP和InDel分别为9 642个（3.09%）和16个（1.52%）；至少在50份里出现的SNP和InDel分别为6 444个（2.06%）和10个（0.95%）；至少在75份里出现的SNP和InDel分别为5 281个（1.69%）和6个（0.57%）（图3）。

基于此，本研究选择至少在50份里出现的SNP和InDel可作为可靠分子标记，它们所对应的基因共有3 742个。

2.4 抗性相关位点与基因分析

本研究利用日本落叶松活动期和休眠期的转录组进行差异基因分析，鉴定到1 7 828个差异表达基因（p＜0.05）（图1）。将这些差异表达基因与带有可靠分子标记的3 742个基因进行韦恩分析（图1），发现有2 569个基因既在活动期和休眠期差异表达，又携带可靠分子标记（图1、4）。

对这2 569个基因进行GO注释：注释到“对真菌的反应（G0：0009620）”，“对真菌的防御反应（GO：0050832）”，“对共生真菌的反应（GO：0009610）”和“对真菌防御反应的调控（GO：1900150）”的基因有26个，其中在活动期上调的有17个（表4）；注释到“对细菌的反应（GO：0009617）”，“对细菌的防御反应（GO：0042742）”和“对细菌防御反应的调控（GO：1900424）”的基因有50个，其中在活动期上调的有31个（表5）；注释到“对卵菌的反应（GO：0002239）”，“对卵菌的防御反应（GO：0002229）”和“对卵菌防御反应的调控（GO：1902288）”的基因有5个（表6）；注释到“对线虫的反应（GO：0009624）”，“对线虫的防御反应（GO：0002215）”，“对昆虫的反应（GO：0009625）”和“对昆虫的防御反应（GO：0002213）”的基因有8个（表7）；注释到“对病毒的反应（GO：0009615）”，“细胞对病毒的反应（GO：0098586）”，“病毒的传播（GO：0046794）”，“对病毒的防御反应（GO：0051607）”，“病毒诱导的基因沉默（GO：0009616）”和“病毒在多细胞宿主中的转运（GO：0046739）”的基因有11个（表8）。

3 讨论

本研究在CDS区共鉴定到515 935个SNP位点。注释后发现，39.45%的SNP位点为同义突变，即这些突变不引起氨基酸改变，这与芒果（Mangifera odorata （Huani））和木瓜（Caricapapaya L．）中同义突变SNP的分布频率一致。同义突变SNP本身不会改变蛋白质序列，因此编码蛋白的稳定性仍可以维持，这可能构成了日本落叶松、芒果和木瓜等植物适应复杂生存环境的遗传基础。

本研究共鉴定到1 056个InDel位点。与SNP相比，InDel发生在CDS区更可能会导致蛋白质结构和功能的改变，从而影响生物体的性状。所以为了维护蛋白质序列和结构的稳定，InDel发生在CDS区上的数量比SNP少，这与在哺乳动物和果蝇（Drosophila melanogaster）中的研究结果一致。

分子标记位点可以揭示林木种群间及种群内的遗传多样性。李培等人利用SRAP标记对来自中国的29个红椿（Toona ciliata Roem.）种源及1个澳大利亚种源进行遗传多样性分析，发现在总的遗传变异中，种源间分化占79.26%，种源内分化仅占20.74%。使用来自不同地区的日本落叶松的转录组数据，本研究发现日本地区的分子标记位点最多（78 219个，41 .42%），说明日本种群的遗传多样性较为丰富。日本是日本落叶松的原生地，随后被引种到不同地区，包括中国、英国。引种后，和原生地存在时空上的隔离，使得不同地区日本落叶松之间的基因交换频率降低。因此，这四个地区共有的位点相对较少，仅有12 879个（6.82%）（图2）。

本研究通过比较活动期与休眠期转录组，发现在3 742个可靠位点所在基因中有2 569个差异表达，其中有101个与真菌、细菌、卵菌、病毒、昆虫和线虫抗性反应相关。这些基因不仅构成了日本落叶松适应生存环境的分子基础，可能也在日本落叶松防御病原物攻击的过程中发挥重要作用；另外，这些基因中存在的SNP和InDel可以作为分子标记在日本落叶松抗性育种中进一步研究和利用。

在真菌感染植物细胞时，植物细胞会通过激活特定的防御机制。本研究共鉴定到17个与真菌相关且在活动期表达量高的基因（表4），其中PUX2存在着错义突变SNP。拟南芥（Arabidopsisthaliana（L．）Heynh．）PUX2的一种突变会降低白粉菌（Golovinomyces orontii）在拟南芥上的繁殖。因此，日本落叶松PUX2的SNP突变也可能影响日本落叶松对真菌的抗性。

本研究共鉴定到31个与细菌相关且在日本落叶松活动期表达量高的基因（表5）。有研究发现，EIN3和EIL1突变后拟南芥增强了对丁香假单胞菌（Pseudomonas syringae）的抵抗力。具有类似情况的基因还有NFXL1，PUB13等。这些结果表明在活动期表达量高的基因在日本落叶松对细菌的抗性反应中可能发挥负调控作用。针对这些负调控抗性的基因，可以通过人工合成miRNA的方式干扰它们的表达以增强日本落叶松的抗性，这种方法在拟南芥、烟草（Nicotianatabacum L．）及番茄（Solanum lycopersicum L．）等植物上的抗病毒研究中已有应用。

本研究鉴定到8个与昆虫和线虫相关且在活动期和休眠期差异表达的基因（表7），其中包括和线虫相关的两个基因ZIF/和MTPc2。已有研究发现，接种根结线虫（Meloidogyne spp.）后拟南芥ZIF/的表达量降低，而MTPc2的表达量在接种后第一周表达先升高，随后降低。这些数据表明，ZIF1和MTPc2在日本落叶松对线虫的防御中发挥一定的作用。

4 结论

通过分析来自不同地区的158份日本落叶松转录组，共获取了515 935个SNP位点和1 056个InDel位点，它们分布在35 827个基因上。根据不同地区转录组中的位点数量，推测日本地区的日本落叶松遗传多样性较为丰富。至少在50份样品中出现的非同义突变SNP位点有6 444个，InDel位点为10个，它们分布在3 742个基因上，可以作为可靠的分子标记进行利用。活动期和休眠期的转录组比较后，发现带有可靠分子标记的2 569基因差异表达；GO注释后发现其中101个基因与植物对真菌、细菌、卵菌、病毒、昆虫和线虫的抗性反应有关。这些结果不仅为日本落叶松全基因组关联分析和全基因组选择育种提供了分子标记，也为利用转基因和基因编辑手段进行遗传育种提供了候选基因。

（责任编辑：张研）

基金项目：科技创新2030 -重大项目（2022ZD0401602; 2022ZD0401705）