Range-wide population genomics of the spongy moth, Lymantria dispar (Erebidae): Implications for biosurveillance, subspecies classification and phylogeography of a destructive moth

AbstractThe spongy moth, Lymantria dispar, is an irruptive forest pest native to Eurasia where its range extends from coast to coast and overspills into northern Africa. Accidentally introduced from Europe in Massachusetts in 1868–69, it is now established in North America where it is considered a highly destructive invasive pest. A fine-scale characterization of its population genetic structure would facilitate identification of source populations for specimens intercepted during ship inspections in North America and would enable mapping of introduction pathways to help prevent future incursions into novel environments. In addition, detailed knowledge of L. dispar’s global population structure would provide new insight into the adequacy of its current subspecies classification system and its phylogeographic history. To address these issues, we generated >2,000 genotyping-by-sequencing-derived SNPs from 1,445 contemporary specimens sampled at 65 locations in 25 countries/3 continents. Using multiple analytical approaches, we identified eight subpopulations that could be further partitioned into 28 groups, achieving unprecedented resolution for this species’ population structure. Although reconciliation between these groupings and the three currently recognized subspecies proved to be challenging, our genetic data confirmed circumscription of the japonica subspecies to Japan. However, the genetic cline observed across continental Eurasia, from L. dispar asiatica in East Asia to L. d. dispar in Western Europe, points to the absence of a sharp geographical boundary (e.g., the Ural Mountains) between these two subspecies, as suggested earlier. Importantly, moths from North America and the Caucasus/Middle East displayed high enough genetic distances from other populations to warrant their consideration as separate subspecies of L. dispar. Finally, in contrast with earlier mtDNA-based investigations that identified the Caucasus as L. dispar’s place of origin, our analyses suggest continental East Asia as its evolutionary cradle, from where it spread to Central Asia and Europe, and to Japan through Korea., MethodsMoth sampling The bulk of spongy moth specimens were collected during the summers of 2017 and 2018, using milk-carton type pheromone-baited traps. Additional samples (whole moths or parts thereof) from regions not fully covered by our network of traps were provided by colleagues who had collected them in the context of independent studies; with the exception of a few samples, collection dates for these were recent (overall, 93% of the moths used were collected between 2013 and 2018, and most specimens were males [only 0.8% of females]) DNA extraction and sequencing For DNA extraction, we sampled one antenna and three legs from each moth. These were frozen in liquid nitrogen and ground using a Retsch MM 200 mixer mill (Retsch technology, Haan, Germany). Then, DNA was extracted with the DNeasy 96 Blood & Tissue Kit (Qiagen, Carlsbad, CA, USA) following the manufacturer's instructions, with the exception of an additional RNase A treatment before the addition of buffer AL/ethanol (4 µL of 100 mg/mL Rnase A; 5 min digestion at room temperature). DNA concentration and purity of the extracts were assessed using a NanoDrop 8000 spectrophotometer (Thermo scientific, Waltham, MA, USA). Samples were diluted to 10 ng/μL prior to library construction. Libraries were prepared based on a genotyping-by-sequencing (GBS) protocol using the restriction enzymes PstI and MspI (Poland et al., 2012). Individuals were barcoded with unique sequences and pooled in multiplexes of 96 individuals per library. Moths from the same sampling site were randomized in the different libraries to reduce the chances of artifactual library effects being interpreted as a biological pattern. Library preparation and sequencing on Ion Torrent Proton P1v2 chips were carried out at the Genomic Analysis Platform of Université Laval, Quebec City, Canada (for a detailed description of the method, see Abed et al. 2019). GBS data The data consists of 1445 sequence individuals from 65 populations and sampled in 25 different counrties, The files Picq_EvolApp_SpongyMothBiosafe_Table_1_20221107.xlsx  and Picq_EvolApp_SpongyMothBiosafe_Figure_population_localisation_20221018.docx give the details for each studied populations (localisation, effective etc.). Raw data 12 plates of 96 individuals; sequencing in 2018; sequencer ouput filename  Cusson_p*_c01.fastq.txt.gz 4 plates of 96 individuals; sequencing in 2019; sequencer ouput filename Cusson_AGM_p*_c01.fastq.gz 1 plate de 96 individuals but only 10 individusals considered in the present project (popultion from Georgia); sequencing in 2020; sequencer ouput filename I.Giguere_BioSAFE_AGM_c01.fastq.gz The file providing the barcodes is named Picq_EvolApp_SpongyMothBiosafe_SampleBarcodes_20221221.xls. Filtered data Picq_EvoApp_SpongyMothBiosafe_datafiltered.vcf For filtering details, please see section 2.4 in the article related to this dataset., Usage notesThe different filtering procedures were carried out using VCFtools v0.1.16 (Danecek et al., 2011), and the resulting VCF file was converted to file formats suitable for each subsequent analysis using PGDSpider v2.1.1.5 (Lischer & Excoffier, 2012). Danecek, P., Auton, A., Abecasis, G., Albers, C.A., Banks, E., DePristo, M.A., Handsaker, R.E., Lunter, G., Marth, G.T., Sherry, S.T., McVean, G., Durbin, R., & 1000 Genomes Project Analysis Group. (2011). The variant call format and VCFtools. Bioinformatics, 27(15), 2156‑2158. DOI: 10.1093/bioinformatics/btr330 Lischer, H.E.L., & Excoffier, L. (2012). PGDSpider : An automated data conversion tool for connecting population genetics and genomics programs. Bioinformatics, 28(2), 298‑299. DOI: 10.1093/bioinformatics/btr642

摘要 舞毒蛾（Lymantria dispar）是一种突发性森林害虫，原产于欧亚大陆，其分布范围横跨整个大陆并蔓延至北非。1868至1869年随人类活动从欧洲意外引入美国马萨诸塞州后，该虫现已在北美定殖，并被视为极具破坏性的入侵害虫。对其种群遗传结构进行高精度解析，将有助于识别北美船舶检疫中截获标本的源种群，并可追溯其入侵路径，助力防范未来向新环境的入侵事件。此外，深入了解舞毒蛾的全球种群结构，可为当前亚种分类系统的合理性及其系统发育历史提供新的研究视角。为解决上述问题，本研究对来自25个国家/3大洲65个采样点的1445份当代标本进行测序，获得了超过2000个基于测序分型（genotyping-by-sequencing, GBS）得到的单核苷酸多态性（Single Nucleotide Polymorphisms, SNPs）位点。通过多种分析方法，我们鉴定出8个亚种群，可进一步划分为28个分组，实现了该物种种群结构解析前所未有的分辨率。尽管将这些分组与当前公认的3个亚种进行对应存在挑战，但遗传数据证实了日本亚种（Lymantria dispar japonica）的分布范围仅限于日本。然而，我们在欧亚大陆全域观察到的遗传梯度，从东亚的舞毒蛾亚洲亚种（L. dispar asiatica）到西欧的舞毒蛾指名亚种（L. d. dispar），表明此前提出的这两个亚种之间存在显著地理界线（如乌拉尔山脉）的观点并不成立。值得注意的是，来自北美和高加索/中东地区的蛾类种群与其他种群间的遗传距离足够显著，足以支持将其划分为舞毒蛾的独立亚种。最后，与此前基于线粒体DNA（mtDNA）的研究将高加索地区定为舞毒蛾起源地不同，本研究的分析结果表明东亚大陆才是其演化起源地，该虫随后从中亚扩散至欧洲，并经朝鲜半岛传入日本。

方法

蛾类采样
本研究的大部分舞毒蛾标本采集于2017和2018年夏季，采用奶盒型性信息素诱捕器进行采集。对于诱捕网络未完全覆盖的区域，我们向曾在独立研究中采集过标本的同行获取了额外样本（完整蛾体或其组织片段）；除少数样本外，这些额外样本的采集时间均较近（总体而言，本研究使用的93%的蛾类标本采集于2013至2018年间，且绝大多数为雄性个体，雌性个体仅占0.8%）。

DNA提取与测序
DNA提取时，我们从每头蛾体上采集1根触角和3条足。将组织样品置于液氮中冷冻后，使用Retsch MM 200混合球磨仪（Retsch technology, 德国哈恩市）进行研磨。随后按照制造商说明书，使用DNeasy 96血液与组织试剂盒（Qiagen，美国加利福尼亚州卡尔斯巴德市）提取DNA，仅在添加AL缓冲液/乙醇前增加了一步RNase A处理步骤（加入4 μL浓度为100 mg/mL的RNase A，室温下消化5分钟）。使用NanoDrop 8000分光光度计（Thermo Scientific，美国马萨诸塞州沃尔瑟姆市）评估提取物的DNA浓度与纯度。将样品稀释至10 ng/μL后进行文库构建。文库制备基于测序分型（genotyping-by-sequencing, GBS）方案，使用限制性内切酶PstI和MspI（Poland等，2012）。每个个体使用独特序列进行条形码标记，每个文库以96个个体为一组进行混合。将同一采样点的蛾类个体随机分配至不同文库，以降低将人工文库效应误判为生物模式的概率。文库制备及Ion Torrent Proton P1v2芯片测序工作均在加拿大魁北克市拉瓦尔大学基因组分析平台完成（详细方法参见Abed等，2019）。

测序分型数据
本数据集包含来自25个国家65个种群的1445份序列个体样本。文件Picq_EvolApp_SpongyMothBiosafe_Table_1_20221107.xlsx 与 Picq_EvolApp_SpongyMothBiosafe_Figure_population_localisation_20221018.docx详细记录了各研究种群的相关信息（采样地点、有效种群大小等）。

原始数据
· 2018年测序的12块96孔板样本，测序输出文件名为Cusson_p*_c01.fastq.txt.gz
· 2019年测序的4块96孔板样本，测序输出文件名为Cusson_AGM_p*_c01.fastq.gz
· 2020年测序的1块96孔板样本（仅10份个体纳入本研究，为来自格鲁吉亚的种群），测序输出文件名为I.Giguere_BioSAFE_AGM_c01.fastq.gz
提供条形码信息的文件名为Picq_EvolApp_SpongyMothBiosafe_SampleBarcodes_20221221.xls。

过滤后数据
Picq_EvoApp_SpongyMothBiosafe_datafiltered.vcf。关于过滤流程的详细信息，请参见本数据集关联论文的第2.4节。

使用说明
本研究使用VCFtools v0.1.16（Danecek等，2011）完成不同的过滤步骤，并使用PGDSpider v2.1.1.5（Lischer & Excoffier，2012）将生成的VCF文件转换为适配各后续分析的文件格式。
参考文献：
Danecek, P., Auton, A., Abecasis, G., Albers, C.A., Banks, E., DePristo, M.A., Handsaker, R.E., Lunter, G., Marth, G.T., Sherry, S.T., McVean, G., Durbin, R., & 1000 Genomes Project Analysis Group. (2011). The variant call format and VCFtools. Bioinformatics, 27(15), 2156‑2158. DOI: 10.1093/bioinformatics/btr330
Lischer, H.E.L., & Excoffier, L. (2012). PGDSpider : An automated data conversion tool for connecting population genetics and genomics programs. Bioinformatics, 28(2), 298‑299. DOI: 10.1093/bioinformatics/btr642