Studying genetic population structure to shed light on the demographic explosion of the rare species Barbitistes vicetinus (Orthoptera, Tettigoniidae).

Insect outbreaks usually involve important ecological and economic consequences for agriculture and forestry. The short-winged bush-cricket Barbitistes vicetinus Galvagni & Fontana, 1993 is a recently described species that was considered rare until ten years ago, when unexpected population outbreaks causing severe defoliations across forests and crops were observed in north-eastern Italy. A genetic approach was used to analyse the origin of outbreak populations. The analysis of two mitochondrial regions (Cytochrome Oxidase I and II and 12S rRNA-Control Region) of 130 samples from the two disjunct ranges (Euganean and Berici Hills) showed high values of haplotype diversity and revealed a high geographical structure among populations of the two ranges. The high genetic variability observed supports the native origin of this species. In addition, results suggest that unexpected outbreaks are not a consequence of a single or few pestiferous haplotypes but rather the source of outbreaks are local populations which have experienced an increase in each area. The recent outbreaks have probably appeared independently of the genetic haplotypes whereas environmental conditions could have affected the outbreak populations. These findings contribute to a growing understanding of the status and evolutionary history of the pest that would be useful for developing and implementing biological control strategies for example by maximizing efforts to locate native natural enemies.

Introduction

Herbivorous insect outbreaks usually involve important ecological and economic consequences for agriculture and forestry. In most cases, it concerns the accidental introduction of an exotic species into a new geographical area due to international trade and human movement [1-3]. In this scenario, the species is not under biological control by natural enemies, leading to a quick spread of the population. Instead, the insect outbreak of an indigenous species is less common but more enigmatic, due to some alterations of the biotic and abiotic components [4].

Until a decade ago, among bush-crickets belonging to the genus Barbitistes, (Orthoptera, Tettigoniidae), occasional outbreaks were reported only for B. constrictus Brunner von Wattenwyl, 1878 in conifer forests of central and eastern Europe, for B. ocskayi (Charpentier, 1850) in broadleaf forests of north-eastern Italy [5-7] and for B. serricauda (Fabricius, 1798) in vineyards of north Italy. In spring 2008, the first outbreak of the species B. vicetinus Galvagni & Fontana, 1993 was recorded in a restricted area of north-eastern Italy (Euganean Hills) [7]. Since then, within two disjunct distribution ranges, the Euganean Hills first and the Berici Hills later, the outbreak areas have progressively enlarged causing heavy damage to forests and neighboring crops (mainly grapes and olives). Severe defoliations have been recorded in the most serious infestations, rising to nearly 90% of canopy loss [8, 9]. Moreover, the outbreaks are also a source of annoyance to people living close to the attacked areas due to the tendency of bush-crickets to invade streets and gardens [7, 10] (S1 Fig).

Interestingly, the bush-cricket was first described just in 1993 as a rare endemic species of north-eastern Italy, being found only in small confined hilly areas [11]. From then on, all the records of B. vicetinus have been reported only in this area. Only since 2005 this species started to be more common in some localities and easy to be found also without bioacoustics technics. The recent discovery of the species (Galvagni & Fontana, 1993) and increase in these north-eastern areas affected by outbreaks during the past two decades has sparked a debate about the possible exotic origin of B. vicetinus. Recently, as a consequence of the outbreaks, some authors have thoroughly studied important aspects of the biology and ecology of B. vicetinus [7, 9, 12, 13] but no data regarding population genetics is so far available. Genetic studies based on the use of high-resolution DNA markers allow to examine the structure of insect populations, identify haplotypes, reconstruct current or past patterns of gene flow and provide information on the origin and expansion routes of the insect [14]. This information could be useful for improving the knowledge concerning the eruptive species, such as colonizing capacity, adaptability, behaviour and demographic history of a population. [15-17]. Several population genetics studies have been conducted with the aim of highlighting and discovering the sources and colonization routes of outbreak populations of alien pests [e.g. 18–23]. However, studies on the genetic structure of native pests have generally received less attention [4, 24].

Among the genetic tools to be exploited, maternally-inherited mtDNA is widely used since it is relatively conserved compared to some nuclear genes, and is thus suitable when searching for historical processes [25]. In this study, we used a combination of two mitochondrial markers, namely Cytochrome Oxidase I and II (COI-tRNALeu-COII) and 12S rRNA-Control Region (12S-CR) to study the genetic structure and demographic history of B. vicetinus populations from the two disjunct ranges where it currently outbreaks (Euganean Hills and Berici Hills, northern Italy). The COI and COII are considered slowly evolving genes of the mitochondrial protein-coding genes [26]. The non-coding control region is in contrast, the most variable segment in the maternally inherited mtDNA [27, 28] and its use in genetic analysis can alleviate biases of coding regions [29]. Both markers have been widely used to conduct population genetic surveys in Orthopteran species including the Tettigoniidae family [30-35].

By studying the genetic structure of the bush-cricket this study aims to answer three main questions:

is there any evidence for an alien origin of the bush-cricket?

are B. vicetinus populations genetically distinct?

are either of the recent outbreaks caused by a single outbreak population followed by a spread pattern or have they derived from multiple local populations?

This information can help understand the factors related to the outbreak events of B. vicetinus and could provide insight to assist management of this outbreak pests.

Materials and methods

Study area

The study area was located in north-eastern Italy (Veneto Region) and included the two disjunct ranges where outbreaks have occurred (Euganean Hills and the Berici Hills) and a third range where no outbreaks have been reported to date and the species is rare (Lessini Mountains). The three ranges are separated from each other by cultivated and inhabited alluvial plains devoid of woody vegetation (Fig 1B).

Fig 1

Parsimony network and geographic distribution of B. vicetinus haplotypes.

A) Haplotypes network based on the combined dataset (COI-tRNALeu-COII and 12S-CR). Each haplotype is represented by a circle, and the area of the circle is proportional to its frequency. Lines within haplotypes circles indicate the proportions shared between collection areas. The colours represent differences in geographic distribution, and small black dots symbolize missing intermediate or unsampled haplotypes. Codes indicate the haplotype ID, reported in Table 1. B) Map showing the proportional geographic distribution of bush-cricket haplotypes across sampled populations. Map republished from [62] under a CC BY license, with permission from [Regione del Veneto–L.R. n. 28/76 –Formazione della Carta Tecnica Regionale], original copyright [2020].

Table 1

Collection sites of B. vicetinus populations analysed and descriptive statistics of each population.

		Site
		ID	Latitude	Longitude	N. samples	Haplotypes	N. haplotypes	H	π (%)
Euganean Hills	North	1	45°20’46.4"N	11°44’4.2"E	1	A6(1)	7	0.64 +/- 0.12	0.04
		2	45°22’30.6"N	11°40’38.7"E	10	A6(7), F6(2), A9(1)
		3	45°21’54.3"N	11°38’52.7"E	8	A1(1), A6(4), A8(2), H6(1)
		4	45°20’59.9"N	11°41’52.8"E	1	G6(1)
	Center	5	45°19’14.8"N	11°45’57.8"E	1	A6(1)	7	0.72 +/- 0.09	0.09
		6	45°17’38.1"N	11°37’22.5"E	2	A6(2)
		7	45°18’14.6"N	11°40’26.8"E	1	D7(1)
		8	45°18’9.07"N	11°40’51.7"E	9	A6(1), A11(1), A12 (2), B6(1), E4(1), E6(3)
		9	45°18’12.4"N	11°43’46.6"E	3	E6(1), A6(2)
		10	45°17’45.6"N	11°39’35.3"E	3	A6(3)
		11	45°18’9.7’’N	11°46’48.5’’E	1	A6(1)
	South	12	45°14’53.8"N	11°44’39.5"E	12	A6(8), A10 (2), C6(2)	4	0.50 +/- 0.13	0.05
		13	45°14’41.1’’N	11°40’6.2’’E	2	A6(1), A13(1)
		14	45°14’55.1"N	11°44’39.3"E	3	A6(3)
Berici Hills	North	15	45°28’20.7’’N	11°35’27.5’’E	16	I1(10), I2(1), I3(1), I15(2), J6(1), L9(1)	9	0.80 +/- 0.06	0.10
		16	45°26’41.2"N	11°34’42.1"E	1	J6(1)
		17	45°27’59.7"N	11°32’52.9"E	8	I2(2), J1(1), J5(1), J6(3), K2(1)
	Center	18	45°25’31.5"N	11°33’02.8"E	1	I5(1)	8	0.90 +/- 0.04	0.13
		19	45°24’00.3"N	11°31’39.7"E	15	A1(1), I1(3), I5(1), I6(3), J1(3), J5(1), J8(2), J16(1)
	South	20	45°25’14.1"N	11°26’59.0"E	17	I5(17)	3	0.28 +/- 0.12	0.02
		21	45°24’54.8’’N	11°28’59.9’’E	3	I1(2), I14(1)
Lessini Mountains		22	45°29’12.3"N	11°02’14.6"E	7	N17(4), O1(3)	3	0.71+/- 0.07	0.24
		23	45°36’14"N	11°25’18"E	5	I1(5)

The Euganean Hills cover an elliptical area of approximately 180 km2 and comprise approximately 100 hills of volcanic origin emerging from the alluvial sedimentary plain, with the highest elevation of approximately 600 m above sea level. They are characterised by numerous narrow and deep valleys, steep hills with different sun exposures and many microclimatic conditions that influence the vegetation. The annual mean temperature is approximately 12°C and precipitation ranges from 700 to 900 mm [36] even if they are extremely variable due to the inner geomorphological variability.

The Berici Hills are an isolated plateau situated on the southern Vicenza plain, they cover almost 160 km2 with a maximum elevation of approximately 450 m above sea level [37]. The climate of the area is characterised by an annual rainfall of 958 mm, and average daily temperature of -1°C in January and 23°C in July [38].

The Lessini Mountains are a triangular-shaped tableland, which occupies some 800 km2, at the transition between the Fore-Alps and the River Po Plain [39]. The mountain group reaches over 2000 m above sea level and is characterised by multiple valleys and long ridges that descend towards the plain.

Insect sampling

Adults of B. vicetinus were collected, using a sweep net, throughout the spring of three successive years (2015, 2016, and 2017). Specimens from the north, centre and south part of both outbreak areas (the Euganean and Berici Hills) as well as from two sites on the Lessini Mountains were collected from bush and tree canopies on marginal side of the forests (Table 1). To avoid sampling relatives, bush-cricket specimens were collected by sweeping an area at least 100 m2 at each sampling site. After capturing, samples were immediately kept in 95% ethanol and taken to the laboratory where they were morphologically identified and stored in individual vials at −20°C until DNA extraction.

DNA extraction, PCR amplification and sequencing

The genomic DNA was extracted from tissue samples taken from the hind femora of each specimen separately, according to a previously described salting-out protocol [40]. Two mitochondrial DNA regions were chosen for amplification: a fragment including the 3’ region of Cytochrome Oxidase I, tRNA-Leu and the 5’ region of Cytochrome Oxidase II (COI-tRNALeu-COII) and a fragment of the 12S rRNA and the Control Region (12S-CR). Two specific primers were designed for amplifying the COI-tRNALeu-COII fragment: Vicetinus_P7F (5’-ACCTGTTCTTGCAGGAGC-3’) and Vicetinus_P4R (5’-TCCACAGATTTCTGAACATTG-3’). The universal primers SR-J14610 (5’-ATAATAGGGTATCTAATCCTAGT-3’) and T1-N18 (5’-CTCTATCAARRTAAYCCTTT-3’) [41] were used to amplify the 12S-CR fragment. Amplifications were performed in 20 μl reactions (1x PCR Go Taq Flexi buffer—Promega, 2.5 mM MgCl2, 0.1 mM dNTPs, 0.5 μM for each primer, 0.5 U of Taq polymerase—Promega, 2 μl DNA template). Thermal cycling conditions for the fragment including the COI-tRNALeu-COII were 5 min at 96°C followed by 35 cycles at 96°C for 1 min, 56°C for 1 min, 72°C for 1:30 min, and a final extension of 72°C for 5 min. For the 12S-CR amplification, the thermal profile followed the conditions described by Eweleit et al. [33] consisting of 2 min at 92°C followed by 35 cycles with a denaturation step of 92°C for 2 min, an annealing step of 52°C for 30 sec, and extension step of 60°C for 3 min, with a final extension of 72°C for 7 min.

PCR products were checked through electrophoresis on 1.0% agarose gels stained with SYBR® (Invitrogen), purified using Exonuclease and Antarctic Phosphatase (GE Healthcare) and sequenced at the BMR Genomics Service (Padua, Italy). Primers used for amplification were also used for sequencing.

Data analysis

DNA sequence chromatograms were quality checked, manually corrected when necessary, and aligned using MEGA X [42]. Low-quality regions found at the beginning and end of each sequence were trimmed, while low-quality sequences were not included in the analysis.

Haplotype and nucleotide diversity, as well as the pairwise genetic distances between populations, were calculated with Arlequin 3.5 [43] using a Kimura 2-parameters model. The presence of population differentiation was also tested by conducting exact tests of population differentiation with 100,000 steps in Markov chain, with 10,000 dememorization steps. To compare the partition of genetic variability among sampled populations, an analysis of molecular variance (AMOVA) [44] was performed using Arlequin.

A haplotype parsimony network with a probability cut-off of 95%, was reconstructed using the TCS 1.21 software [45] and PopART 1.7 [46] and used for depicting the geographical relationships among haplotypes. Ambiguous connections (loops) were resolved using approaches from coalescent theory based on three criteria: frequency, network location and geography [47, 48].

The study of the past demographic history of the species was inferred using Arlequin 3.5 through the Tajima’s D and Fu’s Fs tests [49, 50] and the mismatch distributions of the pairwise genetic differences [51]. Populations at demographic equilibrium or decreasing in size should provide significant positive D and Fs values with a multimodal distribution of pairwise differences, whereas populations that have undergone a sudden demographic expansion usually show significant negative D and Fs values with a unimodal distribution [51, 52]. The sudden expansion model was tested through analysis of the sum of square deviations (SSD) and raggedness index (r) representing the modality of distribution, obtaining the corresponding P values with a parametric bootstrap approach (10,000 replicates).

For those populations that did not deviate from sudden expansion (p > 0.05) the time since expansion was calculated considering the relationship τ = 2*u*t (where τ = age of expansion measured in units of mutational time, u = mutation rate per sequence and per generation, and t = number of generations since the expansion; [53]). The expansion time was then obtained by dividing the estimate of τ by the product of the sequence length in base pairs and the mutation rate per nucleotide (twice the per-lineage substitution rate; u = 2 μk) in percentage per year. Two substitution rates were used: the classical rate of 2.3% divergence per Myr [54] frequently used in Orthoptera [55-59] and the global mtDNA rate of 2.6% divergence per Myr proposed by Papadopoulou et al. [60].

Results

Sampling

A total of 152 individuals were collected from the three natural ranges of B. vicetinus (Euganean Hills, Berici Hills and Lessini Mountains). Unfortunately, only 12 specimens were collected from the Lessini Mountains due to the rarity of the species in this range (Tab. 1).

The two fragments of the mitochondrial genes, COI-tRNALeu-COII and 12S-CR were successfully amplified and sequenced in 130 samples, with an average of 19.6 specimens for each of the three areas within the two main ranges (north, centre, and south of both the Euganean and Berici Hills). After quality assessment and trimming the sequences, high-quality sequences 799 bp long for the COI-tRNALeu-COII and 898 bp long for the 12S-CR fragment were obtained. Sequences of the COI-tRNALeu-COII fragment were translated with Transeq (EMBOSS: http://www.ebi.ac.uk/Tools/emboss/transeq/index.html) to exclude the presence of any nuclear mitochondrial pseudogenes. This bioinformatics tool translates nucleic acid sequences to their corresponding peptide sequences and identifies stop codons. Since these pseudogenes (NUMTs) are characterized by the accumulation of in-frame stop codons and indels [61] the absence of stop codons in the protein sequence can allow to exclude presence of NUMTs.

A total of 13 variable sites, including 8 parsimony informative sites and 5 singleton sites, were identified by the alignment of COI-tRNALeu-COII sequences while the 12S-CR alignment showed 6 parsimony informative sites and 7 singleton sites (13 variable sites in total). Sequences of each haplotype of COI-tRNALeu-COII and 12S-CR obtained in this study are available through GenBank accession numbers MW405351- MW405381.

Partition homogeneity test confirmed that COI-tRNALeu-COII and 12S-CR fragments bear a homogeneous signal (P = 0.28), allowing data to be pooled for further analyses. A final dataset including 130 concatenated sequences of 1697 bp was obtained. Further analyses were conducted considering the combined dataset.

Genetic variability and population structure

The diversity indexes for the concatenated dataset ranged between 0.28 and 0.90 for the haplotype diversity (H) and between 0.02 and 0.24 for the nucleotide diversity (π). Among populations from both the main outbreak areas, those from the north and the centre part of the Berici Hills showed the highest H and π values (H = 0.92, π = 0.13 and H = 0.80, π = 0.10 respectively) while populations from the south part of the Berici Hills showed the lowest variability. The distribution of haplotypes among all populations analysed and other summary statistics are shown in Table 1.

The presence of population differentiation was confirmed by the tests of population differentiation (P>0.001). When B. vicetinus genotypes were grouped based on the six main collection sites (north, centre, and south of both Euganean and Berici Hills), the locus-by-locus AMOVA revealed a significant geographic structure between populations of the two main disjunct ranges (P = 0.0014) (Euganean vs Berici Hills). To avoid bias as a result of the small sampling size of Lessini Mountains, populations from this area were not included in this analysis. The analysis showed that 40% of the variation was explained by differences among groups (Euganean and Berici Hills), whereas about 48% of genetic variation was explained within populations (Table 2).

Table 2

Analysis of molecular variance AMOVA for the combined data set (COI-tRNALeu-COII and 12S-CR).

Structure	Source of variation	Variance %	Fixation indices	P-value
Grouping by geographical range (Euganean and Berici Hills)	Among groups	40.32	FCT = 0.403	0.0014
	Among populations within groups	11.75	FSC = 0.197	<0.001
	Within populations	47.92	FST = 0.520	<0.001

Amova was calculated among populations of Barbitistes vicetinus divided according to the outbreak areas.

TCS Network

TCS Network of the combined dataset revealed the presence of 32 haplotypes of which 15 were exclusive to samples from the Euganean Hills, 13 were exclusive to the Berici Hills and 2 regarded samples only from the Lessini Mountains (Fig 1). Only 2 haplotypes were shared by samples from different ranges: A1 included samples from the Euganean and Berici Hills, and J1 samples from both the Berici Hills and Lessini Mountains. The network evinced a geographical separation among haplotypes, with samples from the Euganean Hills showing a star-like pattern and samples from the Berici Hills connected among them (Fig 1A).

Going deeper into the network’s characteristics, the most common haplotype, A6, included 34 samples coming from all three parts of the Euganean Hills (north, centre and south). Ten rare haplotypes including samples exclusively from the Euganean Hills were connected, separated by only one mutational step, to A6 in a star-shape. Five of them (A9, A11, B6, G6, and H6) were represented by only one sequence, six (C6, A10, F6, A8) by two sequences, and one (E6) by 4 sequences. In addition, haplotype A1, which included samples from both the Euganean and Berici Hills, was also connected to A6 separated from it by only one mutational step.

The second most common haplotype, I1, was composed of 20 sequences:15 belonged to samples from the three geographical populations of Berici Hills (north, centre and south) while 5 were represented by samples from site 23 (Lessini Mountains). This haplotype (I1) was separated from the most dominant haplotype (A6) by two mutational steps. A star-like pattern could also be observed in this part of the network, with haplotype I1 surrounded by nine haplotypes, eight of which (I2, I3, I5, I6, I14, I15, K2 and L9) exclusively represented by samples from the Berici Hills. Among them, haplotype I5 was shared by 19 samples from the southern and central populations of the Berici Hills while haplotypes I2, I3, I15, K2 and L9 were rare (3 sequences each at most) and included only samples for the north part of the Berici Hills. Two haplotypes, O1 and N17, were exclusive to samples from site 22 (Lessini Mountains) the furthest site from the outbreak areas (Fig 1).

The parsimony networks of the single markers (COI-tRNALeu-COII and 12S-CR) showed a similar pattern, with most haplotypes including exclusively samples from only one disjunct distribution range (Euganean Hills or Berici Hills or Lessini Mountains) (S2 and S3 Figs).

Past demographic events

Tajima’s D and Fu’s Fs tests were applied in populations from the Euganean and Berici Hills in order to check for past demographic events. The null hypothesis of neutrality was rejected in populations from the Euganean Hills (D = -2.05; P = 0.02 and Fs = -12.51; P < 0.001), suggesting a past population expansion after a period of low effective sample size (Table 3). Berici populations showed only significant negative Fu’s Fs value (Fs = -6.62; P = 0.004), suggesting that B. vicetinus populations of this hilly area did not conform to the theory of neutral evolution (Table 3).

Table 3

Statistic summary of the past demographic events analysis of Barbitistes vicetinus populations.

Area	Tajima’s D	Fu’s Fs	SSD	r	τ (confidence interval 95%)	Expansion time (ka) 2,3% subst.rate	Expansion time (ka) 2,6% subst.rate
Euganean Hills	-2 .05*	-12 .51*	0 .009	0.069	0.96 (0.34–1.67)	12.4 (4.4–21.6)	10.7 (3.8–18.7)
Berici Hills	-0 .23	-6 .62*	0 .002	0.030	2.14 (0.54–3.62)	27.6 (7–46.7)	24 (6.1–40.3)

Tajima’s (D) and Fu’s neutrality test (Fs) mismatch distribution analysis under a sudden expansion model and time since expansion calculated for mitochondrial populations of the Euganean and Berici Hills considering the combined data base. SSD: sums of squared deviations; r: raggedness index, τ: age of expansion measured in units of mutational time. Expansion time shows as 1,000 years ago (ka).

* Significant at p < 0.05.

The mismatch distribution plots of both ranges were smooth and unimodal, revealing that these populations were undergoing population expansion (S4 Fig). They were characterised by the following observed means: 1.25 and 1.91 for the Euganean and Berici Hills, respectively. Moreover, for both the Berici and Euganean populations, the computed SSD and raggedness index values did not reject a sudden expansion model (S4 Fig), and in particular raggedness values were low (Table 3).

Estimation of the expansion time showed that populations from the Euganean Hills lineage started to expand about 10,700 years ago (with a 2.6% substitution rate) and 12,400 years ago (with a 2.3% substitution rate) (Table 3). Populations from the Berici Hills probably expanded from about 24,000 to 27,600 years ago with 2.6% and 2.3% substitution rates, respectively (Table 3).

Discussion

This study presents the first population genetic analysis of the species B. vicetinus in its outbreak areas and reveals the presence of a high geographical structuring among populations of the two outbreak ranges analysed (Euganean and Berici Hills).

Populations from both ranges showed high values of haplotype diversity, a typical characteristic of ancestral populations [63-65], supporting that B. vicetinus is a native species. Conversely, populations of invasive alien species are traditionally thought to have reduced genetic variation relative to their source populations because of genetic founder effects linked to small population size during the introduction and establishment phases of an invasion [66]. Furthermore, the genetic variability found clustered according to geographical ranges is in contrast with the possibility that B. vicetinus may be an invasive species accidentally introduced.

Even though outbreaks have been reported on both the Euganean and Berici Hills, different haplotypes have been found in these two distribution ranges. Ninety-four percent of haplotypes were exclusive to a single distribution range (i.e. either Euganean or Berici Hills) and only one haplotype, scarcely represented, was found in samples from both hilly areas. These results suggest that outbreaks are not a consequence of a single or few pestiferous haplotypes but rather that the source of outbreaks is due to local populations which experienced a demographic increase in each area. Thus, it seems that outbreaks have appeared independently from the genetic origin, as also found in some studies which indicate that outbreak events are often more affected by environmental conditions than by genetic characteristics of the local population [4, 67–69].

The geographical separation between populations of the Euganean and Berici Hills, observed in the network and confirmed by population differentiation test and AMOVA, indicates a limited gene flow among populations. In phytophagous insects, dispersal capacity, geographical or reproductive barriers, host plant, and habitat fragmentation are reported as the main drivers of genetic structure [20, 70, 71]. The limited dispersal ability of this flightless species could have favoured the lack of gene flow among its distribution ranges. B. vicetinus has been reported to show a low dispersal ability [9], as well as other ground-dispersing species that move only relatively short distances, such as 100–200m during their whole life (e.g. Pholidoptera griseoaptera, [72]). In addition, the mostly lowland areas between the Euganean and Berici Hills, with the presence of agricultural fields and the absence of woody vegetation, might be hostile areas for bush-cricket survival and could have acted as a geographical barrier limiting the effective dispersal of the species. Furthermore, spatial configuration of habitats is another factor influencing genetic structure, that affects mainly species with limited dispersal ability [73]. In B. vicetinus outbreak areas, it has been observed that habitat loss and the presence of patchy areas of the non-host alien tree Robinia pseudoacacia play an important role in reducing population density and dispersion [9]. These factors, coupled to the low mobility of the pest, could have a synergic action that might explain the high level of differentiation observed among populations across its distribution ranges.

High haplotypic diversity, low levels of sequence divergence (nucleotide diversity) and a star-like phylogeny were observed in populations of B. vicetinus from the Euganean Hills and to a lesser extent from the Berici Hills (Fig 1; Table 1). This pattern is consistent with what is expected for populations that have experienced past demographic expansions [52, 74]. These results were highlighted by the neutrality tests and the unimodal mismatch distribution in both hilly areas supporting that, besides the current and ongoing outbreaks, populations from the Euganean and Berici Hills underwent an expansion after a period of low effective sample size (S4 Fig). The bush-cricket population of the Euganean Hills experienced a postglacial expansion starting approximately 10,700–12,400 years ago after the Last Glacial Maximum (LGM; 21,000 years ago) whereas the expansion process in populations on the Berici Hills could have occurred at the end of the LGM. During this glacial period, climate effects were minimal in the Euganean Hills, with thermophilic vegetation serving as a refuge for several species [75, 76]. Accordingly, B. vicetinus populations could have survived during the climatic oscillations, and exploited these hilly areas with potentially suitable environmental conditions throughout the late Pleistocene. Once climatic conditions were favourable, at the end of the last ice age, populations might have experienced an expansion, shaping the present genetic structure of the species. The slow movement of the species and lack of host plants in the lowlands during these periods probably prevented dispersal of the species between the two hilly areas, favouring differentiation of the mitochondrial haplotypes. Thus, both historical (e.g., expansion after LGM) and recent changes (e.g. habitat loss) have contributed to determining the genetic structure and diversity of the bush-cricket in its outbreak range, resulting in genetically structured populations.

Further studies increasing the sampling size, mostly in areas where the species is rare (e.g. Lessini Mountains), combined with analysis of other markers could help obtain a better picture of B. vicetinus population structure.

The current findings contribute to a growing knowledge of the status and evolutionary history of the pest. Here, by removing doubt about the origin of the B. vicetinus, we have achieved an important step in understanding this native species. Shedding light on the origin of a species is also important for its biological control. Given that natural enemies and their host tend to coevolve, biological control programmes often rely on the use of parasitoids present in the native areas. Failure to identify the correct origin of a pest may lead to the use of unsuitable species as biocontrol agents with negative effects on control programmes [20, 77]. Maximizing efforts to locate native natural enemies (native parasitoids), as the egg parasitoid reported in Ortis et al. [13] will be useful for control programmes. Finally, although B. vicetinus is currently in outbreak status, control programme strategies should consider the intrinsic vulnerability of an endemic species. Its low dispersal ability and small and fragmented distribution range could lead this native species to quickly shift from outbreaking to endangered.

Barbitistes vicetinus specimens and damage.

(TIF)

Click here for additional data file.

Parsimony network and geographic distribution of Barbitistes vicetinus haplotypes.

Network obtained with the COI-tRNALeu-COII data set.

(PDF)

Click here for additional data file.

Network obtained with the COI-tRNALeu-COII data set. B) Network obtained with the 12S-CR data set. Network obtained with the 12S-CR data set.

(PDF)

Click here for additional data file.

Mismatch distribution under the population expansion model of Barbitistes vicetinus populations in the disjunct ranges.

(PDF)

Click here for additional data file.

1 Mar 2021

PONE-D-20-40749

Studying genetic population structure to shed light on the demographic explosion of the rare species Barbitistes vicetinus (Orthoptera, Tettigoniidae)

PLOS ONE

Dear Dr. Martinez-Sañudo,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

Please submit your revised manuscript by Apr 10 2021 11:59PM. If you will need more time than this to complete your revisions, please reply to this message or contact the journal office at plosone@plos.org. When you're ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). You should upload this letter as a separate file labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. You should upload this as a separate file labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. You should upload this as a separate file labeled 'Manuscript'.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter. Guidelines for resubmitting your figure files are available below the reviewer comments at the end of this letter.

If applicable, we recommend that you deposit your laboratory protocols in protocols.io to enhance the reproducibility of your results. Protocols.io assigns your protocol its own identifier (DOI) so that it can be cited independently in the future. For instructions see: http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols

We look forward to receiving your revised manuscript.

Kind regards,

Bi-Song Yue, Ph.D

Academic Editor

PLOS ONE

Journal Requirements:

Please review your reference list to ensure that it is complete and correct. If you have cited papers that have been retracted, please include the rationale for doing so in the manuscript text, or remove these references and replace them with relevant current references. Any changes to the reference list should be mentioned in the rebuttal letter that accompanies your revised manuscript. If you need to cite a retracted article, indicate the article’s retracted status in the References list and also include a citation and full reference for the retraction notice.

When submitting your revision, we need you to address these additional requirements.

1. Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at

https://journals.plos.org/plosone/s/file?id=wjVg/PLOSOne_formatting_sample_main_body.pdf and

https://journals.plos.org/plosone/s/file?id=ba62/PLOSOne_formatting_sample_title_authors_affiliations.pdf

2.  Thank you for stating in your Funding Statement:

"This research was partially supported by the project DOR1881327/18 - University of Padua to Luca Mazzon. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript."

Please provide an amended statement that declares *all* the funding or sources of support (whether external or internal to your organization) received during this study, as detailed online in our guide for authors at http://journals.plos.org/plosone/s/submit-now.  Please also include the statement “There was no additional external funding received for this study.” in your updated Funding Statement.

Please include your amended Funding Statement within your cover letter. We will change the online submission form on your behalf.

3. We note that Figure 1 in your submission contain map images which may be copyrighted. All PLOS content is published under the Creative Commons Attribution License (CC BY 4.0), which means that the manuscript, images, and Supporting Information files will be freely available online, and any third party is permitted to access, download, copy, distribute, and use these materials in any way, even commercially, with proper attribution. For these reasons, we cannot publish previously copyrighted maps or satellite images created using proprietary data, such as Google software (Google Maps, Street View, and Earth). For more information, see our copyright guidelines: http://journals.plos.org/plosone/s/licenses-and-copyright.

We require you to either (1) present written permission from the copyright holder to publish these figures specifically under the CC BY 4.0 license, or (2) remove the figures from your submission:

3.1.    You may seek permission from the original copyright holder of Figure 1 to publish the content specifically under the CC BY 4.0 license.

We recommend that you contact the original copyright holder with the Content Permission Form (http://journals.plos.org/plosone/s/file?id=7c09/content-permission-form.pdf) and the following text:

“I request permission for the open-access journal PLOS ONE to publish XXX under the Creative Commons Attribution License (CCAL) CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). Please be aware that this license allows unrestricted use and distribution, even commercially, by third parties. Please reply and provide explicit written permission to publish XXX under a CC BY license and complete the attached form.”

Please upload the completed Content Permission Form or other proof of granted permissions as an "Other" file with your submission.

In the figure caption of the copyrighted figure, please include the following text: “Reprinted from [ref] under a CC BY license, with permission from [name of publisher], original copyright [original copyright year].”

3.2.    If you are unable to obtain permission from the original copyright holder to publish these figures under the CC BY 4.0 license or if the copyright holder’s requirements are incompatible with the CC BY 4.0 license, please either i) remove the figure or ii) supply a replacement figure that complies with the CC BY 4.0 license. Please check copyright information on all replacement figures and update the figure caption with source information. If applicable, please specify in the figure caption text when a figure is similar but not identical to the original image and is therefore for illustrative purposes only.

The following resources for replacing copyrighted map figures may be helpful:

USGS National Map Viewer (public domain): http://viewer.nationalmap.gov/viewer/

The Gateway to Astronaut Photography of Earth (public domain): http://eol.jsc.nasa.gov/sseop/clickmap/

Maps at the CIA (public domain): https://www.cia.gov/library/publications/the-world-factbook/index.html and https://www.cia.gov/library/publications/cia-maps-publications/index.html

NASA Earth Observatory (public domain): http://earthobservatory.nasa.gov/

Landsat: http://landsat.visibleearth.nasa.gov/

USGS EROS (Earth Resources Observatory and Science (EROS) Center) (public domain): http://eros.usgs.gov/#

Natural Earth (public domain): http://www.naturalearthdata.com/

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Yes

**********

2. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: Yes

**********

3. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: No

Reviewer #2: Yes

**********

4. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

5. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: he authors present analyses of mtDNA variation intended to shed light on the origin of a recent pest insect outbreak. The analyses are sound as are the conclusions. I have two main comments and a few minor issues outlined below.

line 169 I think this paragraph should be deleted since the tree is not mentioned (maybe on line 379?) or presented anywhere. Nor is it needed given that assumptions of bifurcating relationships is likely violated (which is presumably the justification for the TCS network).

lines 309-313 One major issue seems to be neglected: how do the estimates of the timing of demographic change compare to the ideas presented in the introduction that the change has occurred very recently? The inferences from molecular data indicate a post-Pleistocene time frame for demographic expansion, but the introduction (line 54) indicates the outbreak has occurred over the last decade. As presented, the attempt to time the demographic events seems to indicate that the authors are thinking of the very recent outbreak, though I suspect, given the nature of the data and the analyses, that only historical events can be detected. I suggest the authors re-arrange the questions and methods to separate these issues and to be more clear on why they are interested in timing the demographic events (which don't have anything to do with the last decade). Temporal sampling with 100's to 1000's of bi-parentally inherited markers might be required to detect changes in Ne over the very near term.

minor comments

line 41 this is not clear - what is meant by "genetic origin"?

line 44 how does mtDNA help with biocontrol strategies? This idea could be expanded by including some of the logic that appears in the discussion.

line 52 I think that some more precise language is needed rather than "ecological balance", a phrase that has had a controversial and tortuous history in ecology.

line 97 DO the authors have the requiste geographical sampling to answer this question of alien origin?? It seems that some basis for comparisons is required for this.

line 114 Change "expositions" to "exposures"

lines 106-121 refer to map fig 1

line 169 I think this paragraph should be deleted since the tree is not presented and not mentioned (maybe on line 379?) or presented anywhere. Nor is it needed given that assumptions of bifurcating relationships is likely violated (which is presumably the justification for the TCS network)?

line 181 delete "a" before coalescent

line 192 delete "occurred"

line 216 it might be useful to describe what you are looking for, i.e. what distinguishes numts from mtDNA genes?

line 229 delete the percent sign

line 379 what aspects are unresolved? The data presented seem to me to provide a fairly clear picture. Further, why "phylogenetic" resolution? The network clearly demonstrates shared haplotypes indicative of gene flow or, more likely, incomplete lineage sorting for which phylogenetic methods are inappropriate - see comment above (given assumption in phylogenetic methods of bifurcating relationships and extinct ancestral states).

FIg 1 top: what do the lines within haplotype circles signify? I think the authors might be indicating proportions shared between regions (haplotypes A1 and I1) or among subregions (i.e. A6, I5) but this is not clear from the colouration.

Reviewer #2: A great application of these techniques. It would be nice to see similar application to other pest species.

I have only minor questions:

(1) How significant is the limited geographical mobility to your interpretation? Would you not come to similar conclusions on species that can disperse but usually show site fidelity?

(2) I would like to see some discussion of application to other species.

(3) What light can such techniques shed on the efficacy of geographically-based pest control methods (such as quarantines)?

**********

6. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: Yes: Shripad Tuljapurkar

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email PLOS at figures@plos.org. Please note that Supporting Information files do not need this step.

25 Mar 2021

Reviewer #1

The authors present analyses of mtDNA variation intended to shed light on the origin of a recent pest insect outbreak. The analyses are sound as are the conclusions. I have two main comments and a few minor issues outlined below.

line 169 I think this paragraph should be deleted since the tree is not mentioned (maybe on line 379?) or presented anywhere. Nor is it needed given that assumptions of bifurcating relationships is likely violated (which is presumably the justification for the TCS network).

Ok. Thank you. We deleted this paragraph.

lines 309-313 One major issue seems to be neglected: how do the estimates of the timing of demographic change compare to the ideas presented in the introduction that the change has occurred very recently? The inferences from molecular data indicate a post-Pleistocene time frame for demographic expansion, but the introduction (line 54) indicates the outbreak has occurred over the last decade. As presented, the attempt to time the demographic events seems to indicate that the authors are thinking of the very recent outbreak, though I suspect, given the nature of the data and the analyses, that only historical events can be detected. I suggest the authors re-arrange the questions and methods to separate these issues and to be more clear on why they are interested in timing the demographic events (which don't have anything to do with the last decade). Temporal sampling with 100's to 1000's of bi-parentally inherited markers might be required to detect changes in Ne over the very near term.

The analyses of past demographic events were conducted to understand how past events could have shaped the genetic structure of B. vicetinus populations. The time since the expansion events revealed by the analysis was calculated to infer the historical period in which the expansion could have taken place, being aware that this expansion is different from the recent outbreak events. However, this analysis is secondary in importance compared to the genetic diversity analysis.

To be clearer and avoiding misunderstandings we delete one sentence in the abstract, modified the third question in the introduction (line 101), as well as some sentences in discussion (lines 366-368), paying more attention to keep separated these issues.

line 41 this is not clear - what is meant by "genetic origin"? We changed this unclear word. Please see line 40

line 44 how does mtDNA help with biocontrol strategies? This idea could be expanded by including some of the logic that appears in the discussion. Thank you for your suggestion. We briefly expanded it. Please see lines 43-44.

line 52 I think that some more precise language is needed rather than "ecological balance", a phrase that has had a controversial and tortuous history in ecology. Done. Please see lines 52.

line 97 DO the authors have the requiste geographical sampling to answer this question of alien origin?? It seems that some basis for comparisons is required for this. We agree with the reviewer. However, the distribution range of B. vicetinus is limited to the small area of north-east Italy where the sampling was performed. The species has not been reported elsewhere, so far.

In addition, we have conducted a preliminary phylogenetic study of the species belonging to the Western-Palearctic genus Barbitistes (through barcoding). Sequences of B. vicetinus were different from sequences of other Barbitistes sp. and all the B. vicetinus specimens grouped in a monophyletic clade (unpublished data).

We added more details in the text. Please see line 67-68.

line 114 Change "expositions" to "exposures” Done. Please line 115.

lines 106-121 refer to map fig 1 Done. Please line 111.

line 169 I think this paragraph should be deleted since the tree is not presented and not mentioned (maybe on line 379?) or presented anywhere. Nor is it needed given that assumptions of bifurcating relationships is likely violated (which is presumably the justification for the TCS network)? Done.

line 181 delete "a" before coalescent Done. Please see line 179.

line 192 delete "occurred" Thank you. Done. See line 190.

line 216 it might be useful to describe what you are looking for, i.e. what distinguishes numts from mtDNA genes? Ok. Thank you. We added more information. See line 214-218.

line 229 delete the percent sign Done. Please see line 232.

line 379 what aspects are unresolved? The data presented seem to me to provide a fairly clear picture. Further, why "phylogenetic" resolution? The network clearly demonstrates shared haplotypes indicative of gene flow or, more likely, incomplete lineage sorting for which phylogenetic methods are inappropriate - see comment above (given assumption in phylogenetic methods of bifurcating relationships and extinct ancestral states). We agree with the reviewer. We modified it. Please see line 384-385.

FIg 1 top: what do the lines within haplotype circles signify? I think the authors might be indicating proportions shared between regions (haplotypes A1 and I1) or among subregions (i.e. A6, I5) but this is not clear from the colouration. Thank you for your suggestion

We added this information in the figure legend trying to make it clearer:

Lines within haplotypes circles indicate the proportions shared between collection areas

(Please see lines 263-264)

Reviewer #2

A great application of these techniques. It would be nice to see similar application to other pest species.

I have only minor questions:

(1) How significant is the limited geographical mobility to your interpretation? Barbitistes vicetinus tends to stay close to its hatching site during the whole life cycle, likely due to its low dispersal ability. This characteristic, coupled to the low density of the species until 2008, probably helped to maintain the strong genetic structure of the species. However, in correspondence with the outbreaks, a tendency of the insect to disperse more than before was reported (Cavaletto et al. 2018). During the outbreaks, B. vicetinus have colonized new habitats and moved far from the hatching areas, perhaps to avoid intraspecific competition.

Besides the low mobility of the bush-cricket other factors such as geographical barriers and spatial configuration could have favoured the absence of gene flow between B. vicetinus populations. The lowland areas between the Euganean and Berici Hills, with several agricultural fields and without woody vegetation could have acted as a geographical barrier, notwithstanding the short geographical distance. We have also observed that patchy areas of the non-host alien tree Robinia pseudoacacia have played an important role in reducing population density and dispersion of the pest (Cavaletto et al. 2019).

(1) b) Would you not come to similar conclusions on species that can disperse but usually show site fidelity? The low dispersal ability is undoubtedly a contributing factor to maintaining a high levels of genetic differentiation. Other species that have the ability to move but show fidelity to their life site could also show a strong genetic structure, although it may be less marked as with low-mobility species. Species showing breeding site fidelity could move if necessary (e.g. habitat loss, climate changes) and it could modify the genetic structure.

(2) I would like to see some discussion of application to other species. Similar studies can be applied to a wide range of insect species. Population genetic strategies have been widely used for examining patterns and magnitude of dispersal of several insects over geographic and temporal scales. Regarding invasive species, for example, population genetics studies can be conducted for highlighting and discovering the sources and colonization routes of these alien insects.

The information obtained with genetic analyses can be also useful for an effective conservation of endangered species. Indeed, analysing the genetic diversity patterns of endangered species has become an integral component of many management strategies.

Some examples of the application of these techniques to other species can be found in lines 74-83 and 93-95.

(3) What light can such techniques shed on the efficacy of geographically-based pest control methods (such as quarantines)? The genetic techniques used in this study have allowed us to realize that population outbreaks are not always a consequence of a single or few pestiferous haplotypes. As reported here, any individual of this species, regardless of its genetic haplotype, had an equal chance of producing outbreaks. It follows that control strategies must always be applied (e.g.: strategies based on geographic control).

For species characterized by a reduced movement behaviour throughout their lifespan (such as the bush-cricket here studied), gene flow appears low enough to maintain a high genetic structure among populations. In this case, control programs based on isolated treatments (i.e. at a small geographic scale) could be also effective.

Thus, these results based on genetic techniques strengthen the need to apply and enhance control methods such as quarantines or other strategies based on geographic control.

Submitted filename: Response to Reviewers.docx

Click here for additional data file.

8 Apr 2021

Studying genetic population structure to shed light on the demographic explosion of the rare species Barbitistes vicetinus (Orthoptera, Tettigoniidae)

PONE-D-20-40749R1

Dear Dr. Martinez-Sañudo,

We’re pleased to inform you that your manuscript has been judged scientifically suitable for publication and will be formally accepted for publication once it meets all outstanding technical requirements.

Within one week, you’ll receive an e-mail detailing the required amendments. When these have been addressed, you’ll receive a formal acceptance letter and your manuscript will be scheduled for publication.

An invoice for payment will follow shortly after the formal acceptance. To ensure an efficient process, please log into Editorial Manager at http://www.editorialmanager.com/pone/, click the 'Update My Information' link at the top of the page, and double check that your user information is up-to-date. If you have any billing related questions, please contact our Author Billing department directly at authorbilling@plos.org.

If your institution or institutions have a press office, please notify them about your upcoming paper to help maximize its impact. If they’ll be preparing press materials, please inform our press team as soon as possible -- no later than 48 hours after receiving the formal acceptance. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information, please contact onepress@plos.org.

Kind regards,

Bi-Song Yue, Ph.D

Academic Editor

PLOS ONE

22 Apr 2021

PONE-D-20-40749R1

Studying genetic population structure to shed light on the demographic explosion of the rare species Barbitistes vicetinus (Orthoptera, Tettigoniidae)

Dear Dr. Martinez-Sañudo:

I'm pleased to inform you that your manuscript has been deemed suitable for publication in PLOS ONE. Congratulations! Your manuscript is now with our production department.

If your institution or institutions have a press office, please let them know about your upcoming paper now to help maximize its impact. If they'll be preparing press materials, please inform our press team within the next 48 hours. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information please contact onepress@plos.org.

If we can help with anything else, please email us at plosone@plos.org.

Thank you for submitting your work to PLOS ONE and supporting open access.

Kind regards,

PLOS ONE Editorial Office Staff

on behalf of

Dr. Bi-Song Yue

Academic Editor

PLOS ONE