Invasion genetics of the Asian hornet Vespa velutina nigrithorax in Southern Europe

In 2004, Vespa velutina was first seen in France. Since then, this fierce honey bee predator spread across many countries, giving rise to one of the most phenomenal insect invasions in Europe. An early study in France showed a genetically depauperate population, originating from a single multi-mated queen introduced from China. Here, we further unveil V. velutina invasion genetics in Europe by surveying the Iberian and Italian peninsulas using cytonuclear markers. Our results show that the French population acted as the colonists’ source in Spain, Portugal and Italy, leading to rejecting the hypothesis of multiple introductions from the native range. While Spain and Italy were colonized predominantly by leading-edge expansions from the French core population, in Portugal the invasion started from long-distance jump. Both processes were accompanied by a significant reduction in genetic diversity, with stronger losses for Portugal (Ar = 17.4%; uHe = 42.3%) than for Spain (Ar = 9.0%; uHe = 20.6%) or Italy (Ar = 16.3%; uHe = 26.8%). Signatures of differentiation and population structure, associated to the founding event in Portugal, enabled detection of secondary contact between the front derived from the primary propagule introduced in France and the front derived from the secondary propagule introduced in Portugal. Detection of first-generation migrants in the three countries suggests continuous gene flow that is bringing in new alleles, and this effect is stronger in Portugal, as reflected by a 20.3% increase in allelic richness. Overall, this study provides further insights into the invasion genetics of V. velutina in Europe, which can aid developing strategies to manage this major threat to beekeeping.

expansions from the French core population, in Portugal the invasion started from long-distance jump. Both processes were accompanied by a significant reduction in genetic diversity, with stronger losses for Portugal (Ar = 17.4%; uHe = 42.3%) than for Spain (Ar = 9.0%; uHe = 20.6%) or Italy (Ar = 16.3%; uHe = 26.8%). Signatures of differentiation and population structure, associated to the founding event in Portugal, enabled detection of secondary contact between the front derived from the primary propagule introduced in France and the front derived from the secondary propagule introduced in Portugal. Detection of first-generation migrants in the three countries suggests continuous gene flow that is bringing in new alleles, and this effect is stronger in Portugal, as reflected by a 20.3% increase in allelic richness. Overall, this study provides further insights into the invasion genetics of V. velutina in Europe, which can Supplementary Information The online version contains supplementary material available at https://doi.org/10.1007/ s10530-022-02730-9.

Introduction
In the past century, numerous insects have successfully colonized new territories around the globe (reviewed in Beggs et al. 2011;Elder and Bell 1998;Kenis et al. 2009;Paine and Millar 2002). Yet, of special concern are social Hymenoptera, with six species in the top-100 IUCN list of the worst alien invaders (http://go.nature.com/qa9z1g). The Asian hornet, Vespa velutina nigrithorax de Buysson (hereafter referred to as V. velutina), has not yet earned a place in the IUCN infamous list. However, the phenomenal spread across western Europe, adversely affecting biodiversity, ecosystem services and even human health, led to its inclusion in the European Union's (EU) black-list of invasive alien species (EU regulation No. 1143. V. velutina shares with many other social Hymenoptera a number of unique biological and life history traits that enhance the probability of successful establishment and subsequent spread in new territories (reviewed in Moller 1996). Polyandry is one such trait and V. velutina is one of the rare polyandrous species in Vespidae (Arca et al. 2015). Depending on the mating frequency in the original habitat, a few (or even a single) founder queens may carry in their spermatheca a representative sample of the genetic diversity of the source population (Ding et al. 2017;Mikheyev et al. 2009;Schmid-Hempel et al. 2007). This trait greatly aids the early stage of propagule establishment by lessening the effects of genetic drift and increasing the likelihood of fast adaptive response to the new selective regimes in the recipient environment (Estoup et al. 2016). Once established, a suit of other traits provides additional adaptive benefits, facilitating spread across the novel habitats, including (i) exceptional dispersal abilities, (ii) high reproductive efficiency, (iii) broad diets and habitat ranges, (iv) effective predator defences, and (v) ability of buffering against environmental changes (reviewed in Beggs et al. 2011;Moller 1996).
While over 30 species of Vespidae were successfully introduced around the world in the last decades (Beggs et al. 2011), Europe was free of these invaders until arrival of a single mated queen of V. velutina (Arca et al. 2015) to the French district of Lot-et-Garonne in 2004 (Arca et al. 2015;Haxaire et al. 2006;Villemant et al. 2006). V. velutina has since expanded across most of France and made its way into nine neighbouring countries. It took six years to the first incursion of V. velutina outside of France. In 2010, the first specimen was captured in the Spanish locality of Amaiur, 180 km south of the introduction point in Nérac, Lot-et-Garonne (López et al. 2011). However, soon after Spain, several other countries registered the entry of V. velutina in their territories, including Belgium  and Portugal (Grosso-Silva and Maia 2012) in 2011, Italy in 2012 , Germany in 2014 (Witt 2015), Great Britain ) and the Netherlands in 2016 (Smit et al. 2018), and Switzerland in 2017 (Poidatz et al. 2018). The recent detection in Hamburg, Germany, far from the leading edge in the Netherlands, represents the northernmost point in Europe (Husemann et al. 2020). While the geographical limits in Europe are yet to be reached (Barbet-Massin et al. 2013, 2018Bessa et al. 2016;Keeling et al. 2017;Robinet et al. 2017), currently, V. velutina occupies a large continuous territory in the western countries, from Portugal to the Netherlands and across the English Channel in Great Britain, and it is spreading eastwards into Italy (see the dynamic map in http://frelonasiatique.mnhn.fr/).
The spread of V. velutina in Europe is causing economic and ecological losses as well as human fatalities (Monceau et al. 2014). While the yearly costs of nest destruction (the main control measure) have exceeded €29 million in the United Kingdom, France and Italy (Barbet-Massin et al. 2020), the damage to wild ecosystems, agriculture and beekeeping is likely much larger but yet to be fully appraised. V. velutina is disrupting invaded ecosystems by competing with the native entomofauna (Monceau et al. 2015) and, primarily, by interfering with pollination of wild plants and crops as it feeds upon many pollinator insects (Monceau et al. 2013a), most notably the western honey bee (Apis mellifera). In Autumn, when the hornet's nest size and protein intake by the sexual brood reach their maxima (Monceau et al. 2013b), hundreds of hornets hover at the front of beehives disrupting foraging, a critical activity at this time of the year as colonies are preparing for winter. This disturbance can lead to depletion of food stores and to depopulation of beehives, enhancing the likelihood of winter colony mortality (Requier et al. 2019).
Spain, Portugal, and Italy are home to 5.4 million colonies, representing 29.4% (as in 2019) of the honey bee stock and 22.3% (as in 2018) of the honey production in the European Union. In a scenario of growing colony losses (Neumann and Carreck 2010), sustaining these figures implies additional expenses with apiary management and thereby a heavier financial burden to the beekeeping industry. The drivers of colony losses are complex and interacting and include a diverse repertoire of alien pests, such as the ectoparasitic mite Varroa destructor, and its vectored viruses, the beetle Aethina tumida, the microsporidia Nosema ceranae, and now V. velutina (Laurino et al. 2020;Martín-Hernández et al. 2018;Mutinelli et al. 2014;Steinhauer et al. 2018). Unfortunately, the early attempts to halt establishment and spread of V. velutina in France were insufficient and inefficient. Also inefficient were the Pyrenees and the Alps in deterring its natural range expansion into the Iberian and Italian peninsulas (Bertolino et al. 2016;Goldarazena et al. 2015;López et al. 2011). On the other hand, the arrival of V. velutina to Viana do Castelo (Minho, Portugal), in 2011, and to A Mariña (Galiza, Spain), in 2012, over 620 and 480 km apart, respectively, from the expanding front in Navarra and Basque Country (Grosso-Silva and Maia 2012; López et al. 2011), likely resulted from a human-mediated jump within Europe, although the hypothesis of independent introductions from the native range cannot be overlooked.
V. velutina has established itself successfully in Iberia, where it occupies a continuous area from Catalonia (Mediterranean) to Lisbon (Atlantic), and it is currently expanding in the Italian peninsula along the western coast (Bertolino et al. 2016;Porporato et al. 2014). The mean expansion rate ranges from 18.3 in Italy (Bertolino et al. 2016) to 37.7 km year -1 in Portugal (Carvalho et al. 2020), suggesting that V. velutina is spreading at a slower pace in these territories than in France, where an invasion speed of up to 82 km year -1 was documented (Robinet et al. 2017). Despite the surveillance and control actions enforced by Reg. (EU) No. 1143/2014 to the EU member states, populations have rapidly grown to large numbers, as reflected by nest densities of 3.5 nests km -2 in Italy (Bertolino et al. 2016) and 5.4 nests km -2 in Portugal (Carvalho et al. 2020).
While there is a wealth of studies on different aspects of the biology, ecology and management of V. velutina in Europe (reviewed in Laurino et al. 2020;Monceau et al. 2014;Turchi and Derijard 2018), there are only a few addressing the genetic aspects of the invasion (Arca et al. 2015;Budge et al. 2017;Granato et al. 2019;Jones et al. 2020). In these studies, mitochondrial and nuclear DNA data are consistent with a founder event starting from a single mated queen (Arca et al. 2015), whose descendants reached Italy by natural expansion (Bertolino et al. 2016;Granato et al. 2019) and Great Britain by humanmediated jump ). Here, we expanded current knowledge on the invasion genetics of V. velutina in Europe by reporting on the results of a mitochondrial (cytochrome C oxidase I) and nuclear (microsatellites) DNA screening of the invading populations in the Iberian and Italian peninsulas. Our aim was (i) to test whether the invading populations derive from the introduced propagule in France or from the native range in Southeast Asia, (ii) to determine the extent of genetic diversity impoverishment relative to the French and Asian populations, and (iii) to infer genetic structure.

Materials and methods
Sampling and DNA extraction A total of 246 adult V. velutina females were collected in Spain (N = 45), Portugal (N = 190) and Italy (N = 11) from founding nests, secondary nests, and from traps placed in apiaries between March and December of 2016 and 2017 (see sampling details in Table S1). Apiaries were selected to be at least 3 km far apart to avoid over-representing the queen genotype by sampling multiple individuals from the same colony. Samples were placed in absolute ethanol and stored at -20°C until molecular analysis. Total DNA was extracted from the thorax of the 246 individuals, each representing a single colony, using the Ron's Tissue DNA Mini Kit (ÒBioron), according to the protocol provided by the manufacturer.
Mitochondrial DNA sequencing and data analysis A sub-sample of 13 individuals from Spain, 22 from Portugal, and 3 from Italy (Table S1), selected in each country to cover a wide geographical range, was amplified at the mitochondrial DNA (mtDNA) 642 bp barcode region of the Cytochrome C Oxidase Subunit I (COI) gene, using the universal primers LCO-1490 and HCO-2198 (Folmer et al. 1994). Polymerase chain reaction (PCR) amplification was carried out to a final volume of 25 lL which contained 0.2 mM of dNTPs (Promega), 2 mM of MgCl 2 , 1X Taq DNA polymerase buffer, 0.2 lM of each primer, 0.625U of GoTaqÒ-Flexi DNA polymerase (Promega) and 10 ng of DNA. PCR was performed in a T100 Thermal Cycler (BioRad TM ) using an initial denaturation step of 948C for 5 min, followed by 35 cycles of 30 s at 94°C, 45 s at 50°C and 1 min at 72°C, and a final extension of 10 min at 72°C. PCR products were sent to STABVIDA Inc. (Portugal) for direct Sanger sequencing in both directions. DNA sequences were checked manually for base calling and aligned using MEGA 7.0.26 (Kumar et al. 2016).
The 38 sequences were aligned manually with 23 sequences of V. velutina downloaded from GenBank using MEGA. To infer the phylogenetic relationships among the 61 samples, the maximum parsimony method was implemented in MEGA. The tree was obtained using the subtree-pruning-regrafting (SPR) algorithm (Nei and Kumar 2000) with search level zero, in which the initial trees were obtained by the random addition of sequences (10 replicates). The resulting tree was rooted with Vespa bicolor (accession number KT257112), Vespa vivax (KT257116) and Vespa affinis (KJ147242).

Microsatellite genotyping and data analysis
The 246 individuals were genotyped for 21 microsatellite loci using the primers developed by others (Arca et al. 2012;Daly et al. 2002;Hasegawa and Takahashi 2002). The forward primers were labelled with 6-FAM, HEX, Atto 550 or Atto 565 fluorescent dyes and the PCR amplifications were performed in three multiplexes (Table S2). The PCR reactions contained 1 lL of DNA (10 ng/lL), 5 lL of Multiplex PCR kit (Qiagen), and 1 lL of a 2 lM primer mix, in a 10 lL final volume. Reactions were performed in a T100 Thermal Cycler (BioRad TM ) using a temperature profile consisting of an initial denaturation step at 95°C for 15 min followed by 40 cycles of 95°C for 30 s, 55°C for 1 min and 72°C for 1 min, and a final extension of 60°C for 60 min. Fragment analysis was performed on an ABI PRISM 3730 9 l automated Genetic Analyzer (Applied Biosystems, Foster, California, USA), using GeneScan 500LIZÒ as the internal size standard, at STABVIDA Inc. (Portugal). The lengths of the fragments were determined using GENEMAPPER 3.7 (Applied Biosystems) and each allele was checked manually. To have an estimate of the amplification error rate, 30 randomly selected samples were genotyped twice and analysed with GIMLET 1.3.3 (Valière 2002). Further scoring errors due to null alleles and allelic dropout were assessed in MICRO-CHECKER v.2.2.3 (Van Oosterhout et al. 2004).
The microsatellite dataset generated for the 246 individuals was merged with an existing dataset containing genotypes of 162 individuals collected in the introduced (France = 83) and native ranges (Chinese provinces of Yunnan = 20 and Zhejiang/ Jiangsu = 30, Indonesia = 21, and Vietnam = 8) of V. velutina (Arca et al. 2015). To enable dataset merging, allele scores of each microsatellite locus were harmonized between laboratories by genotyping 10 DNA samples used by Arca et al. (2015).
Allele frequencies, observed number of alleles (Na), effective number of alleles (Ne), observed heterozygosity (Ho), expected heterozygosity (He) and unbiased expected heterozygosity (uHe) were computed per locus and country using GENALEX 6.5 (Peakall and Smouse 2012) whereas allelic richness (Ar) was computed using the package HIERFSTAT (Goudet 2005) in RStudio v1.2.5033 (RStudio Team, 2015). The diversity metrics uHe, calculated using the standard approach of Nei (1978), and Ar, calculated using the rarefaction method (El Mousadik and Petit 1996), are adjusted for unequal sample size. Differences in Na, Ar, and uHe between population samples were assessed by the Wilcoxon's signed-rank test, also implemented in RStudio. Genetic differentiation among populations from the invaded and native ranges was estimated by pairwise F ST values computed in ARLEQUIN 3.5.2.2 (Excoffier et al. 2005), using 10 000 random permutations to assess significance. Whenever needed, statistical significance levels were adjusted for multiple comparisons using the sequential Bonferroni procedure (Rice 1989).
To identify the most likely source of invaders in Spain, Portugal and Italy, the 246 individuals were assigned to the Asian and French reference populations from Arca et al. (2015) using the Bayesian-based assignment test (Rannala and Mountain 1997) and the Monte-Carlo resampling algorithm of Paetkau et al. (2004), with 10 000 replicates, in GENECLASS 2 (Piry et al. 2004). The assignment test estimates the likelihood that an individual belongs to a population and assigns the population with the highest likelihood as the origin of the sampled individual. Considering that not all potential source populations have been sampled, first-generation migrants in each of the three invaded countries were identified using L_home as the statistical criterion for likelihood estimation (Piry et al. 2004). A threshold of 0.01 was set for type I error.
To assess genetic relationships among individuals and whether there is genetic structure in Europe, a model-free principal coordinate analysis (PCoA) was run in GENALEX 6.5 (Peakall and Smouse 2012), using the option Genetic Distance (Codom-Genotypic) for codominant data, and a model-based Bayesian clustering analysis was run in STRUCTURE 2.3.3 (Pritchard et al. 2000) using the merged dataset (N = 408). The STRUCTURE analysis was performed using the admixture ancestry and correlated allele frequency models with the unsupervised option to estimate the membership proportion (Q-value) to each ancestral cluster (K) for each individual. The program was set up for 750 000 Monte Carlo Markov chain iterations, after an initial burn-in of 250 000, and run for 20 independent replicates for each K (from 1 to 5). Q-plots were post-processed online with CLUM-PAK (Kopelman et al. 2015).

Mitochondrial diversity
Mitochondrial DNA (mtDNA) variation was analysed using sequence data of the 642 bp COI barcode region obtained from 38 individuals from Spain, Portugal, and Italy. A single haplotype, named by Arca et al. (2015) as F, was identified in these individuals, as shown in the maximum parsimony phylogenetic tree (Fig. 1). The European samples were placed in the same clade, together with two samples (JQ780450 and JQ780451) from Zhejiang, and with a close phylogenetic proximity to a sequence from Yunnan (JQ780457), further supporting the claim of eastern China as the origin of V. velutina propagules in Europe (Arca et al. 2015).

Microsatellite diversity
A total of 21 microsatellite loci were screened in the 246 individuals from Spain, Portugal, and Italy. Five of these loci were initially excluded from further analysis due to inconsistent amplification. The remaining 16 loci were successfully amplified with virtually no missing data (see genotypes in Table S1 and missing data in Table S2). According to the tests run in GIMLET and MICRO-CHECKER, there were barely any signs of allelic dropout (Tables S3 and S4). In contrast, null alleles were detected for all loci, but VMA8, and in all countries, but Italy (Table S4). Yet, the 15 loci did not overlap within and between the European and the Asian populations, suggesting that the loci themselves are not faulty. Instead, it is possible that the null alleles are an artifact from joining samples (collected throughout wide geographical areas and in different years) by country, causing a Wahlund effect. Support for this hypothesis comes from the geographically-narrow sample collected in Italy and from the sub-sample collected in 2016 in Viana do Castelo (the place of V. velutina introduction in Portugal) with no null alleles (Table S4). Accordingly, we did not discard any locus and performed the ensuing analysis using the 16 loci.
The mean (± SE) number of alleles (Na, Ne), allelic richness (Ar), and heterozygosities (Ho, uHe) computed for each country of the invaded and native ranges in Europe and Asia are shown in Table 1 (but  see Table S5 for allele frequencies and metrics per locus and population). As expected for a founder event, the genetic diversity of V. velutina was dramatically reduced following its introduction in Europe. The 16 microsatellite loci had a total of 190 different alleles in the eight population samples analysed, of which 186 were present in Asia and 65 in Europe. While the number of private alleles varied between five (Vietnam) and 20 (Zhejiang/Jiangsu) in the native range, in Europe there was only one private allele detected in the Spanish population. In Asia, all loci were polymorphic whereas in Europe two loci were monomorphic (R1-77 in Spain, Portugal, and Italy; R1-80 only in Portugal). Significantly fewer (P-value B 0.0036, Wilcoxon's signed-rank test) alleles were detected in Italy (2.062 ± 0.170), Portugal (2.312 ± 0.198) and Spain (2.438 ± 0.182) than in France (3.938 ± 0.452) or in the Chinese province of Zhejiang/Jiangsu (7.062 ± 1.039).

Source of invaders
The most likely origin of V. velutina invasion in Spain, Portugal, and Italy was identified by analysing all pairwise F ST values and mean probability assignments using the Asian and French populations as potential sources (Table 2). A high level of genetic differentiation was observed for the intercontinental comparisons (mean ± SE; 0.362 ± 0.031), being the lowest F ST values obtained with the Chinese province of b Fig. 1 Maximum Parsimony (MP) phylogenetic tree for the 38 V. velutina individuals from Spain, Portugal, and Italy sequenced in this study. The MP tree was generated from 64 COI sequences with a length of 598 bp (256 nucleotide positions were eliminated because they contained gaps or missing data). The MP was obtained with a tree length of 196. MP inferred a consistency index of 0.844, a retention index of 0.912, and a composite index of 0.770 for parsimonyinformative sites. The tree is drawn to scale, with branch lengths calculated using the average pathway method (Nei and Kumar 2000) and are in the units of the number of changes over the whole 598-bp COI sequence. The analysis included 64 additional sequences from individuals of the invaded range in Europe, South Korea, and Japan, and from the native range in China (Yunnan, Jiangsu/Zhejiang), Vietnam and Indonesia The mean assignment probabilities for the three invasive populations are shown in Table 2. Every single individual sampled in Spain, Portugal, and Italy was assigned to France (Table S6), with a score of 100% and assignment probabilities varying between 0.277 and 0.999 for Spain (0.735 ± 0.031), 0.262 and 0.979 for Portugal (0.695 ± 0.013), and 0.410 and 0.961 for Italy (0.668 ± 0.058). These metrics are considerably higher than those obtained for Zhejiang/ Jiangsu, the second most likely source of origin, with a mean score of 0% and mean assignment probabilities varying between 0.329 ± 0.026 for Spain and 0.434 ± 0.010 for Portugal (Tables 2 and S6). Altogether, the results from F ST and assignment tests suggest that France acted as the most likely source of the colonists that invaded both the Iberian and Italian peninsulas.

Population structure
Genetic relationships of all individuals from Europe and Asia were first investigated using the model-free principal coordinate analysis (PCoA). Individuals from Portugal and Asia fell at the two extremes of the spectrum defined by coordinate 1 axis, which  explains 24.2% of the total variation (Fig. 3). Most individuals from Portugal (96.3%) and Spain (53.3%) formed a tight cluster, which is placed in the threedimensional space further apart from the Asian cluster than from the cluster that joins all individuals from France and Italy (Fig. 3). A pattern of genetic heterogeneity in the invaded range in Europe was also captured by the Bayesian model-based approach implemented in STRUCTURE (Fig. 4). At K = 2, there was a clear separation between invasive and native populations. Further partitioning of variation (K = 3) unveils the existence of genetic structure in Europe, with the individuals from France (mean Q-value ± SE; 0.958 ± 0.007; Table S7) and Italy (0.975 ± 0.008) placed together in one cluster, coloured by orange in the graphical representation of structure (Fig. 4). In contrast, individuals from Portugal were placed together in the blue cluster (0.958 ± 0.013), except for six individuals that exhibited a clear ancestry in the orange cluster (Qvalue C 0.934). Individuals from Spain displayed mean intermediate Q-values (0.528 ± 0.071), with 19 individuals assigned to the orange cluster (0.927 B Q-value B 0.995) and 23 to the blue cluster (0.910 B Q-value B 0.996). The six Portuguese individuals assigned to the orange cluster were also detected by GENECLASS as first-generation migrants (P-value B 0.0019; Table S8) with a putative origin in Spain (7.936 B -log(L) B 10.547). GENECLASS also identified two migrants in Spain (P-value B 0.0006; 9.998 B -log(L) B 11.425) and one in Italy (P-value B 0.0000; -log(L) B 8.965), but they were both originated from the French population (Fig. 4). Interestingly, when the first-generation migrants were removed from the population samples, there was a decrease in diversity in Spain, Portugal, and Italy, which was more pronounced for allelic than for heterozygosity measures, with Ar varying between 20.3% (Portugal) and 13.2% (Italy) and uHe varying between 3.3% (Portugal) and 2.0% (Spain. See Fig. 5 and Table S9).

Discussion
The mitochondrial and nuclear markers concurrently employed herein suggest that the invasions of V. velutina into the Iberian and Italian peninsulas stemmed from the particularly successful population of France, rather than from multiple independent introductions from the native range. Convincing evidence comes from the unique mitochondrial  Table S7 for Q-values) for each individual of the invading and native ranges (X-axis). Q-values C 0.90 into the French cluster and the blue dots with Q-values C 0.85 into the Portuguese cluster. The arrows point to first-generation migrants, identified using L-home as the statistical criterion for likelihood estimation, with a threshold of 0.01 for type I error (see Table S8  and Italy (IT), when firstgeneration migrants are removed from each population sample. Loss uHe = 1 -(uHe d / uHe r ) and Loss Ar = 1 -(Ar d /Ar r ), being allelic richness (Ar) and unbiased heterozygosity (uHe) calculated across loci for each sample with firstgeneration migrants discarded from the dataset (d) and retained in the dataset (r) haplotype identified in the area under study, congruent with other mitochondrial surveys undertaken in Europe (Arca et al. 2015;Budge et al. 2017;Granato et al. 2019;Jones et al. 2020;Takahashi et al. 2019), and with the hypothesis of a single mated-queen propagule introduced into Europe ( Arca et al. 2015).
The nuclear data further supports the assertion of a French source, as microsatellite alleles identified in the peninsular populations are a subset of the alleles found in France, which in turn are a subset of the alleles found in Asia (Arca et al. 2015). A similar finding has been reported by the other single microsatellite surveys conducted up to now outside of France Jones et al. 2020). Furthermore, the F ST values indicate that the peninsular populations are more closely related to the French population than to any other Asian population. This is in agreement with the assignment test, which assigned every single individual to the French population with a 100% score and a high probability. Observational data suggest that V. velutina crossed the French borders, into Spain and Italy, at the northwestern edges of the Pyrenes and the Alps (Bertolino et al. 2016;Goldarazena et al. 2015;López et al. 2011). The STRUCTURE analysis performed herein supports these dispersal pathways, with the individuals from Italy and the Basque Country (Spain) clustering with the French population. In contrast, the founding populations early detected in Portugal, far from the invading front in the Basque Country (Goldarazena et al. 2015;López et al. 2011), could have originated from independent propagules, consistent with mounting evidence that multiple introductions from native sources are common in biological invasions (reviewed in Bossdorf et al. 2005;Dlugosch and Parker 2008;Roman and Darling 2007;Wares et al. 2005). However, our data does not support a scenario of multiple introductions pointing instead to a long-distance jump into Portugal, as has been documented for the United Kingdom Jones et al. 2020). Whether the propagule originated from the core population in France or from the leading edge, at the time located in the Spanish Basque Country, cannot be determined by this dataset. However, there are rumours that the invading propagule hitchhiked with a pine wood cargo imported from the Landes region (France) into Viana do Castelo (Portugal) in 2010 (M. Maia, pers. comm.). This mechanism of long-distance jump has greatly accelerated the invasion of Europe, earlier by expanding the southwesternmost edge in Portugal (Grosso-Silva and Maia 2012) and more recently the northwesternmost edge in the United Kingdom Jones et al. 2020).
Descendants of the Portuguese founding population rapidly established and spread throughout the Atlantic coast southwards and northwards into the neighbouring Spanish province of Galiza (Carvalho et al. 2020). The first sightings of V. velutina in Galiza, only one year after its official arrival in Portugal, occurred in A Mariña, over 200 km north from the founding nucleus in Viana do Castelo and 400 km west from the invading front in the Spanish Basque Country (Goldarazena et al. 2015;López et al. 2011). Considering an average dispersal rate of * 40 km year -1 in the Iberian Atlantic coast (Carvalho et al. 2020), the introduction in A Mariña was most likely assisted by humans. Interestingly, STRUCTURE places the individuals sampled in the region of A Mariña in the Portuguese blue cluster, contrasting with all neighbouring individuals which belong to the French orange cluster. The genetic pattern shown in Fig. 4 suggests that Galiza is a place of secondary contact between two fronts: one stemming from the founding population in Portugal and the other from the natural range expansion along the Atlantic coast in northern Iberia (see the dynamic map in http://frelonasiatique.mnhn. fr/). Whether the individuals sampled in the region of A Mariña with ancestry in the blue cluster derive from the founding propagule detected in Galiza in 2012 or represent recent human-mediated migration with origin in the Portuguese population is uncertain and deserves closer scrutiny.
Founder or range expansion events are usually accompanied by losses of genetic diversity, which can be more or less severe depending on propagule size (Nei et al. 1975) or on the influx of migrants from the core into edge populations (Swaegers et al. 2013). The impact of the founder event on neutral diversity was of intermediate intensity in the French population of V. velutina, as it retaining 85.2% of the allelic richness and 79.9% of the heterozygosity present in the putative source population of China (Arca et al. 2015). This level of diversity loss is in the range documented for numerous animal and plant invasions, with varying propagule pressures, and was claimed to be linked to their invasion success (reviewed in Bossdorf et al. 2005;Dlugosch and Parker 2008;Roman and Darling 2007;Wares et al. 2005). Empirical support for this claim comes from a meta-analysis on plants and animals that found a significant positive effect of genetic diversity on measures of colonization (Forsman 2014).
Although founder diversity is not the only factor determining the fate of introduced propagules, it is a critical one as it will influence the adaptive response of fitness-related traits to the novel environment (Drake and Lodge 2006;Lockwood et al. 2005). Remarkably, owing to the polyandrous mating system in V. velutina, the single multi-mated queen that gave rise to the invasion in Europe was able to bring in an important fraction of neutral, and perhaps fitnessrelated, variation (Arca et al. 2015). This trait, together with other extraordinary life history traits common to social Hymenoptera (reviewed in Beggs et al. 2011;Moller 1996), probably helped establishment and spread of her descendants in Europe, similar to invasions of the little fire ant Wasmannia auropunctata, the eastern honey bee Apis cerana, or the bumble bee Bombus terrestris, which also started from a few introduced queens (Ding et al. 2017;Eloff et al. 2020;Mikheyev et al. 2009;Schmid-Hempel et al. 2007).
More severe was the bottleneck associated to the secondary founder event in Portugal, as reflected by the significant losses of 29.6% and 53.9% in allelic richness and heterozygosity, respectively. Notably, when the values were computed relative to the French source, diversity loss was more pronounced in this secondary founder event (Ar = 17.4%; uHe = 42.3%) than in the original introduction in France (Ar = 14.8%; uHe = 20.1%). Nevertheless, the founding propagule was still able to establish and spread throughout the Atlantic coast of Iberia, although at a slower pace than in France (Carvalho et al. 2020). Given suitable climatic conditions and prey abundance in the area (Carvalho et al. 2020;Villemant et al. 2011), it is possible that this diversity impoverishment might have contributed to deceleration of the invasion.
While microsatellite variation is not necessarily a proxy for fitness-related diversity (Reed and Frankham 2001), it is possible that the sex locus suffered further depletion during the secondary founding event, hampering population growth. Similar to other haplodiploid hymenopteran species, in V. velutina, sex is determined by a multiallelic locus in which females develop from heterozygous genotypes whereas males from hemizygous (Heimpel and De Boer 2008). When allelic diversity is low at the sex locus, there is a higher chance that the queen will produce homozygous eggs, which will develop into diploid males. The problem is that these males represent a heavy burden to the colony because, in addition to usually being sterile, they are produced at the expense of female workers but do not contribute to colony tasks. Studies of nests collected in the United Kingdom and France documented a high frequency of diploid males, which is consistent with low diversity at the sex locus and inbreeding among the members of the small founding population Darrouzet et al. 2015;Jones et al. 2020).
While the invasion in Portugal arose from a secondary founder event, the greater levels of genetic diversity in Spain and Italy, as well as the clustering patterns ( Figs. 3 and 4), suggest a process of natural expansion from the core population in France. Leading-edge populations usually lose diversity at lower rates than founding populations because they are buffered against the effects of genetic drift by gene flow from the core population, especially in cases of organisms with high dispersal rate as for V. velutina (Shi and Chen 2012;Swaegers et al. 2013). Contrary to expectations from ecological niche modelling, the geographical extent of the invasion in Italy is very modest as compared to that of Iberia (Bertolino et al. 2016;Carvalho et al. 2020;Rome et al. 2011;Villemant et al. 2011). While several interacting biotic and abiotic factors (e.g., clime, prey abundance, interspecific competition, and control measures) may explain differential invasion dynamics between the two peninsulas, genetic factors associated to the secondary introduction in Portugal have likely played an important role in population growth in Iberia. In less than 6 years, the primary front expanded throughout the north Atlantic coast in Spain all the way to Galiza where it met the secondary front originating from Portugal (http://frelonasiatique.mnhn.fr/). Detection of six first-generation migrants in Portugal suggests that the primary front is expanding southward and this movement is bringing in new genetic variation, as reflected by the 20.33% increase in allelic richness (Fig. 5).
Regardless of the many details that remain elusive, our genetic analysis suggests that V. velutina invaded north-western Spain and Italy by a process of natural range expansion whereas the invasive process in Portugal stemmed from a long-distance jump. While severe losses of genetic variability from serial bottlenecks did not avert the successful spread of V. velutina in Portugal, the recent secondary contact between the two invading fronts in Galiza is helping restoring diversity levels, which may increase population fitness and boost invasiveness. However, the findings presented herein come from selectively neutral markers and whether these can be used as surrogates for adaptively important variation is a controversial matter (Reed and Frankham 2001). A recent sharp rise in the number of detected nests in Portugal denotes population growth (Carvalho et al. 2020). Whether this observation is due to improved fitness related with the influx of genetic variation carried by the primary expanding wave or to mechanisms involving phenotypic plasticity, epigenetic modifications, de novo mutations, purging of deleterious mutations, and/or increase of additive variance at ecologically important traits remains to be known (see the review of Estoup et al. 2016).
In this study, we unveiled the patterns of neutral diversity in populations of the Italian and Iberian peninsulas, thereby expanding current knowledge on the invasion genetics of V. velutina. However, a more complete understanding of the molecular mechanisms underpinning the phenomenal invasion of this hornet in Europe requires that functionally important diversity is scrutinized and this can be done using highthroughput sequencing. Unfortunately, the whole genome of V. velutina has not yet been sequenced and therefore the genomic resources for this type of study are still lacking.

Conclusion
This study builds upon the seminal analysis of the founding propagule in France (Arca et al. 2015) to greatly expand current understanding on the population genetics aspects underlying the phenomenal invasion of V. velutina in Europe. Invading events in southern Europe have been achieved by leading edge expansion and long-distance jump. Remarkably, the genetic impoverishment exhibited by the invading populations, which was more severe in the secondary founding population in Portugal than in Spain and Italy, did not preclude successful population establishment and spread, encouraging further research into the mechanisms underlying invasiveness in social Hymenoptera. Equally remarkable is the substantial gain in genetic diversity introduced by a low number of first-generation migrants, suggesting that the signature of genetic drift can be rapidly erased by gene flow and this can occur rapidly in species like V. velutina, with high dispersal abilities. Our findings further stress the importance of surveillance and control measures to halt gene flow, and call for more stringent policies within the European Union preventing circulation of propagules, which may lead to secondary founder events favouring establishment and expansion of alien organisms.