Research Article
Print
Research Article
Hidden on mountaintops: phylogeography and evolution of Parnassiana (Orthoptera: Tettigoniidae) in the sky-islands of the southern Balkans
expand article infoNefeli Kotitsa‡, Simeon B. Borissov‡, Dragan P. Chobanov‡
‡ Institute of Biodiversity and Ecosystem Research, Bulgarian Academy of Sciences, Sofia, Bulgaria

Corresponding author: Dragan P. Chobanov ( dchobanov@gmail.com )

Academic editor: Lara-Sophie Dey
© 2026 Nefeli Kotitsa, Simeon B. Borissov, Dragan P. Chobanov.
This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Citation:
ZooBank: urn:lsid:zoobank.org:pub:4C70C13D-C03C-4639-A7D3-12900606F853
Open Access

Abstract

Parnassiana is a micropterous, flightless genus of bush-crickets that inhabits the mountains of the southwestern Balkan Peninsula, displaying a sky-island distribution. It consists of 13 valid species and many populations with unknown taxonomic status. This study explores the phylogeography and evolution of Parnassiana by reconstructing its evolutionary relationships using multilocus DNA data (NAD2, COI and ITS), by estimating the divergence times of its lineages and relating them to significant geological and climatic events, and by delimiting the boundaries between Parnassiana taxa. Phylogenies are compared with morphological characters and bioacoustics data to inform species delimitation conclusions and evolutionary mechanisms. We conclude that the evolution of Parnassiana has been primarily shaped by tectonic and climatic events of the Pliocene (separation of Peloponnese from Central Greece, establishment of Mediterranean climate) and Early Pleistocene (warm and long interglacial periods), which led to allopatric speciation. Secondarily, the Middle and Late Pleistocene were characterized by dispersal, intra-species diversification, and possible gene exchange, adding complexity to the geographic pattern of the genus and causing discrepancies between the gene trees and phenotypic grouping. The evolutionary history of Parnassiana is reflected in its geographic distribution patterns, as isolated mountains host well-defined older lineages, while lineages from large massifs show a more uniform phenotype and complex phylogenetic relationships.

Key words

bush-crickets, interglacial refugia, morphology, mitochondrial phylogeny, Pliocene, Pleistocene, Platycleidini, morpho-acoustic evolution

1. Introduction

Parnassiana Zeuner, 1941 (Orthoptera: Platycleidini) is a genus of micropterous, flightless bush-crickets (Fig. 1) that is found in the mountains of the southwestern Balkan Peninsula. The members of the genus are restricted to the high altitudes (usually above 1500 m) of the isolated mountain summits of the Pindos range in mainland Greece, the mountains of Peloponnese, and the island of Evvoia (Willemse et al. 2018). The northernmost limits of their distribution include the mountains Smolikas and Tymphi (Papigko), and the southernmost locality is Mt. Taygetos in Peloponnese (Willemse et al. 2018). Along this latitudinal span of around 400 km, populations of the genus occur on at least 33 mountain summits (Willemse et al. 2018), where they seem to be isolated from each other and have given rise to many endemic taxa.

Figure 1. 

Parnassiana in their habitat. Numbers, when present, refer to males (1) and females (2). A P. parnon; B P. menalon; C P. chelmos unicolor – Kyllini Mt.; D P. chelmos chelmos, male – Chelmos Mt.; E P. chelmos unicolor, female – Panachaiko Mt.; F P. chelmos deplanata, male – Erymanthos Mt.; G P. fusca; H P. dirphys; I P. gionica; J P. coracis, male – Vardousia Mt.; K P. coracis – Oxia Mt.; L P. tymphrestos – Tymphrestos Mt.; M Parnassiana sp. 1 – Kaliakouda Mt.; N P. panaetolikon; O Parnassiana sp. 3 – Agrafa Mt.; P Parnassiana sp. 5 – Voutsikaki Mt.; Q P. cf. tymphiensis – Smolikas Mt.; R P. tenuis – Tzoumerka Mt.; S P. tymphiensis – Koziakas Mt.; T Parnassiana sp. 2 – Karava Mt.; U Parnassiana sp. 4 – Avgo Mt.

The disjunct and restricted distribution of Parnassiana, coupled with the large number of endemic species found in such a small area, suggests that the mountains hosting the genus are acting as sky islands. Sky islands are isolated mountain systems separated by lowland areas that function as true islands from an evolutionary and biogeographical perspective, providing conditions that promote isolation, diversification, and speciation (Dodge 1943).

Currently, Parnassiana consists of 13 described species and many populations with unclear taxonomic status, due to subtle morphological differentiation and unknown acoustic communication (Willemse and Willemse 2008; Willemse et al. 2018). Willemse and Willemse (2008) note that populations on summits that are part of larger mountain ranges or groups, such as the mountains of South Pindos in Central Greece (e.g. Tymphrestos, Oiti, Kaliakouda, Vardousia), or the Central and North Pindos range (Delidimi, Avgo, Karava, Tzoumerka, Tymphi and more), are especially challenging to delimit due to their subtle morphological differences.

Four species and three subspecies have been described from the mountains of Peloponnese, namely P. fusca (Brunner von Wattenwyl) (Mt. Taygetos), P. parnon Willemse (Mt. Parnon), P. menalon Willemse (Mt. Menalo), and P. chelmos Zeuner (P. ch. chelmos – Mt. Chelmos; P. ch. deplanata (Willemse) – Mt. Erymanthos; P. ch. unicolor (Willemse) – Mts Kyllini and Panachaikon). From Central Greece, seven species are known: P. tymphrestos Zeuner (Mt. Tymphrestos and Oiti, and tentatively Mts Kaliakouda and Helidona), P. coracis (Ramme) (Mt. Vardousia, and tentatively Mt. Oxia), P. gionica La Greca & Messina (Mt. Giona), P. nigromarginata (Willemse & Willemse) (Mt. Akarnanika), P. dirphys Willemse (Mt. Dirphy), P. parnassica (Ramme) (Mts Parnassos and Elikonas), and P. panaetolikon Willemse (Mt. Panaetolikon) (Willemse et al. 2018; Stefanidis et al. 2025). The taxonomic status of some of these species is unclear, especially considering P. coracis, P. tymphrestos and P. panaetolikon, which show significant morphological similarities and have neighboring distributional ranges (Willemse and Willemse 2008). Lastly, two species are known from the North Pindos range – P. tenuis (Heller & Willemse) from Tzoumerka and adjacent mountains (e.g., Chatzi, Lakmos and Kakarditsa), and P. tymphiensis (Willemse) from the mountains of Tymphi, Smolikas, Vasilitsa and Mavrovouni (Willemse et al. 2018). In the central and south parts of Pindos some mountains such as Tzoumerka and East Agrafa have been reported to host a “panaetolikon-like” taxon, while enigmatic populations which could represent new species are reported from many additional mountain peaks, such as Mt. Karava, Avgo, Voutsikaki/Kazarma, Valtou (and its northern peak Gavrogo), and Delidimi (Heller 2006; Willemse and Willemse 2008; Willemse et al. 2018).

Currently, all Parnassiana species are assessed with Threatened categories in the IUCN Red List due to their extremely restricted distribution ranges and declining populations and habitat quality (IUCN 2025). As a result, Parnassiana has recently drawn the attention of ecologists, especially in light of climate change and its effects on alpine ecosystems. A series of papers by Stefanidis et al. (2024, 2025) highlight the ecological and microhabitat preferences of four Parnassiana taxa found in the mountains of Central Greece (P. parnassica, P. tymphrestos, P. gionica, P. coracis). The authors stress on the importance of phylogenetic studies and the clarification of the systematics of the genus, as any changes in the taxonomic status of the taxa would directly affect their conservation status and, consequently, the suggested conservation strategies (Stefanidis et al. 2025).

Phylogenetic analyses based on DNA sequence data represent a powerful tool for drawing the boundaries of species and identifying Evolutionary Significant Units. They are especially helpful for cryptic and morphologically challenging populations (Zhang et al. 2013; Luo et al. 2018), as well as for reconstructing evolutionary histories, and have been applied successfully for numerous Orthopterans (e.g. Song et al. 2015; Mugleston et al. 2018; Borissov et al. 2023). Alone, however, they are not enough for accurately delimiting the taxa of a group. Evidence based on ecological, morphological and behavioral traits among others is crucial for providing diagnostic traits for taxa and for interpreting clades and relationships that derive from molecular phylogeny (Padial et al. 2010; Zhang et al. 2013). Parnassiana is a prime example of a case that would benefit greatly from an integrative molecular, morphological and bioacoustics treatment, due to its subtle morphological and acoustic differentiation, sky-island distribution, and conservation importance (Heller 2006; Willemse and Willemse 2008; Stefanidis et al. 2025).

The aim of this study is to explore the phylogeography and evolution of a highly threatened genus with sky-island distribution in the southern Balkan Peninsula, using a complete dataset covering the whole distribution of Parnassiana, by: 1) reconstructing the internal phylogeny of Parnassiana populations using one nuclear and two mitochondrial markers, 2) estimating the divergence times of the Parnassiana lineages and relating lineage splits to significant geological and climatic events, and 3) delimiting the boundaries between Parnassiana taxa using a combination of molecular phylogenetic methods, morphology and bioacoustics.

2. Material and Methods

2.1. Sampling

Collecting trips were performed in the mountains of Peloponnese, Central Greece, Evvoia, and the Pindos range in the years 2021, 2022 and 2024, leading to a total of 28 mountain summits visited (Fig. 3). Localities were selected based on the list of mountains from which Parnassiana specimens had previously been collected (Willemse et al. 2018), as well as mountain summits lacking records but exhibiting conditions similar to those of known localities, mainly elevations above 1900 m. Maps were created in QGIS Desktop v.3.28.12 (QGIS.org 2023) using an elevation gradient raster from the SRTM database (Consortium for Spatial Information (CGIAR-CSI) 2004, https://srtm.csi.cgiar.org/srtmdata/), and freely available vector data (IUCN Basedata 2024). Additional outgroup material from taxa of Platycleidini, with a focus on putative phylogenetically close lineages, was collected from Albania, Bulgaria, Kazakhstan, North Macedonia, or downloaded from GenBank (www.ncbi.nlm.nih.gov/genbank).

Figure 2. 

Phylogenetic tree of the genus Parnassiana (left panel), inferred from the concatenated NAD2, COI, and ITS sequences using Bayesian Inference analysis. The colors indicate the different phylogenetic lineages and putative taxa. Numbers after the population names and localities correspond to the numbers in Fig. 3. Black dots on nodes indicate Bayesian posterior probability >0.95, while numbers indicate pp<0.95. Lineages corresponding to species delineation entities are marked with bars of different grey scale on the right of the tree, as follows: (A) GMYC; (B) ABGD; (C) ASAP; (D) consensus subjective delimitation based on combined results from this study. In the right panel, male titillators (left) and cerci (right) are visualized for each population.

Figure 3. 

Geospatial distribution of the phylogenetic lineages of Parnassiana. Different major geomorphological units are shown in background colors: North and Central Pindos range – blue; mountains of Central Greece – pink; mountains of the Peloponnese – brown. Branch colors indicate four well supported major groups: Peloponnese clade – brown; South Pindos subclade (SPS) – red; North-Central Pindos subclade (NCP) – blue; Dirphys-Parnassos subclade (DPS) – ochre. The Taygetos clade and Akarnanika subclade with unresolved positions are marked with purple dashed line. The numbers correspond to taxa and populations included in the publication and to the numbers in Fig. 2.

Details for each specimen used in the present study (coordinates, collection data, locality, GenBank accession numbers, etc.) can be found in File S1. Newly collected specimens are deposited in the Institute of Biodiversity and Ecosystem Research, Bulgarian Academy of Sciences (IBER-BAS) and the Natural History Museum of Crete (NHMC). For the morphological study, a small number of dried specimens from the Naturalis Biodiversity Center’s collection were also examined, including type series and topotypes of several species (P. tymphrestos, P. nigromarginata, P. tymphiensis).

2.2. Molecular phylogeny and time estimations

2.2.1. DNA extraction

Total genomic DNA was extracted from hind femora muscles of 69 specimens of Parnassiana, which originate from 28 mountain summits in Greece, using the Invitrogen PureLink Genomic DNA Mini Kit (Thermo Fisher Scientific 2024). DNA was also extracted from own material of 24 representatives of tribus Platycleidini (following Cigliano et al. 2026), and one sequence of Anterastes babadaghi Uvarov was added from GenBank, which were used as outgroups. There was a focus on Platycleidini genera found in the Balkan Peninsula that have morphological and ecological similarities to Parnassiana (Massa and Fontana 2011; Çiplak et al. 2015). The DNA extraction protocol followed the manufacturer’s instructions.

2.2.2. DNA amplification

One nuclear (the internal transcribed spacers 1 and 2, together with the 5.8 S ribosomal RNA between them, ITS1–5.8S–ITS2, henceforth mentioned as ITS) and two mitochondrial (NADH dehydrogenase subunit 2—NAD2 and cytochrome oxidase subunit 1—COI) DNA markers were used. ITS was amplified with the primers WeekF TAGAGGAAGTAAAAGTCG (forward) and WeekR GCTTAAATTCAGCGG (reverse), resulting in one ITS1–5.8S–ITS2 fragment (Weekers et al. 2001). NAD2 was amplified with the primer pair TM-J210 AATTAAGCTAATGGGTTCATACCC (forward) and TW-N1284 AYAGCTTTGAARGYTATTAGTTT (reverse) (Simon et al. 2006). In cases where TM-J210 did not yield satisfactory results, the primer TI-J34 GCCTGATTAAAGGRTTAYYTTGATA was used as forward instead (Simon et al. 2006). The COI fragment was amplified with the primer pair C1-J-1718 GGAGGATTCGGAAATTGATTAGTACC (forward) and TL-2-N-3014 TCTAATGCATTAATCTGCCATCTTA (reverse) (Simon et al. 1994, modified for locusts).

Polymerase chain reactions were carried out in 25 μl volume using the Thermo Fischer Scientific ‘DreamTaq Hot Start Master Mix’ according to the manufacturer’s instructions. Temperature cycling for mitochondrial fragments followed Chobanov et al. (2017), with adaptations for hot-start PCR and slight adjustments. NΑD2 was amplified with initial step at 95°C for 5 min, followed by 35 cycles of denaturation (95 °C for 50 s), annealing (51–59°C for 40 s), and elongation (72°C for 80 s), with a final elongation step at 72°C for 15 min. For the COI fragment, an initial step at 95°C was held for 5 min., followed by 35 cycles, including denaturation (94°C for 40 s), annealing (50°C for 40 s), and elongation (70°C for 1:30 s). Final elongation step was performed at 72°C for 15 min. For the ITS fragment, the protocol by Ullrich et al. (2010) was applied. Purification of PCR products and Sanger sequencing from both 5’ and 3’ ends were performed by Macrogen Europe (Macrogen, Inc., Amsterdam, the Netherlands).

2.2.3. Phylogenetic analyses

Chromatograms were processed, trimmed and assembled using CodonCode Aligner v.8.0.2 (CodonCode, Dedham, MA, USA). Sequence alignments were performed in MEGA X v.11.0.13 (Kumar et al. 2018) using the MUSCLE algorithm. The absence of significant saturation and of stop codons (for the protein-coding sequences) were confirmed in DAMBE v.7.3.32 (Xia 2018).

Nucleotide substitution models were calculated with ModelFinder (Kalyaanamoorthy et al. 2017), using partitions by coding positions (Chernomor et al. 2016), on the IQ-TREE web server (Trifinopoulos et al. 2016; http://iqtree.cibiv.univie.ac.at). Best models were selected under the corrected Akaike Information Criterion (AICc) (File S2).

Phylogenetic analyses were performed on matrices of each genetic marker as well as on the concatenated matrix of all markers, using Maximum likelihood (ML) and Bayesian Inference (BI). The ML analysis was performed in IQ-TREE (Nguyen et al. 2015), applying the calculated models. Bootstrap support was obtained through ultrafast bootstrap with 1000 replicates (Hoang et al. 2018). The BI analysis was performed with MrBayes v. 3.2.7 (Ronquist et al. 2012) using the closest available approximations of the substitution models suggested by ModelFinder (File S2). Parameters for the Bayesian Inference analyses include four simulations of Markov chains and 2 × 106 generations sampling each 100th tree. Stationary distribution of the MCMC parameters was checked with Tracer v. 1.7.1 (Rambaut et al. 2018). The first 25% of trees were excluded as burnin. Resulting trees were visualized in FigTree v.1.4.4 (http://tree.bio.ed.ac.uk/software/figtree).

2.2.4. Estimation of divergence times

Dating was performed with BEAST and BEAUtie v.2.7 (Bouckaert et al. 2019). BEAST was run with datasets of 31 ingroup (Parnassiana) sequences, corresponding to a single sequence per mountain summit per taxon, and 14 outgroups. The dataset, which included the concatenated NAD2–COI–ITS matrix, was split in five partitions to allow differences in clock rates: two for each protein-coding marker (codons 1 + 2 and codon 3), and one for ITS. Two calibration schemes were applied independently, as follows:

a) Biogeographical calibration point. The early stages of the Corinthian Rift development, which led to the opening of the Corinthian Gulf and therefore the isolation of Peloponnese from Central Greece 4–3.5 mya (Rohais and Moretti 2017; Gawthorpe et al. 2018; Hemelsdaël et al. 2021; Wicker et al. 2024), were used as a calibration point. By 3.5 mya a large part of northern Peloponnese was submerged, and the remaining landmass was effectively an island (Fassoulas 2018). This event has been proven to present a barrier promoting allopatric diversification for multiple organisms, such as beetles, land snails, spiders, lizards, and the Orthopteran genus Poecilimon (Kotsakiozi et al. 2012; Gkontas et al. 2016; Kornilios et al. 2016; Psonis et al. 2018; Borissov et al. 2020). Since Parnassiana is a mountainous genus that cannot tolerate lowland climatic conditions, we hypothesize that the initial stages of rifting in that area caused the permanent isolation of the Peloponnesian lineages (De Baets et al. 2016). Hence, we set the BEAST prior using a normal distribution with a mean of 3.8 mya and standard deviation of 0.15 for the age of split between the most recent common ancestor (MRCA) of the mainland clade and their closest Peloponnesian relative based on our phylogenetic tree.

b) Evolutionary clock rate. Two different clock rates for NAD2 were tested: the evolutionary rate of micropterous orthopterans (Chang et al. 2020) (0.014256 subs/s/million years) and the evolutionary rate for Pholidopterini, calibrated based on the mid-Aegean Trench by Çiplak et al. (2022) (0.018 subs/s/million years). These rates were both applied to the codons 1 + 2 partition of NAD2. The clock rates of the other four partitions of the concatenated matrix (including the codon 3 partition of NAD2) remained unlinked and were estimated by BEAST.

A gamma site model with four category counts and a GTR substitution model were applied for all partitions. The Yule speciation process was set as prior. An MCMC chain length of 108 generations was run, sampling each 1000th tree. A lognormal relaxed clock was applied (Drummond et al. 2006), and the site and clock models remained unlinked across partitions, while trees were linked. Effective sample size (ESS) and stationarity were assessed for all parameters with Tracer v.1.7.1 (Rambaut et al. 2018). Results were summarized with TreeAnnotator v.2.7 (Bouckaert et al. 2019) discarding 10% of the trees as burn-in, and the maximum clade credibility tree was visualized with FigTree v.1.4.4 (Rambaut 2006). The position of P. fusca as a sister lineage of the ‘Mainland clade’ was forced as a prior, in order to match the topology supported by the BI phylogenetic tree.

2.2.5. Sequence-based species delimitation

To delimit the taxa, the phylogenetic species concept, as defined by De Queiroz (2007), Fujisawa and Barraclough (2013), and Zhang et al. (2013), was applied through the GMYC, ABGD and ASAP species delimitation methods.

For the GMYC method, built for single-locus data (Fujisawa and Barraclough 2013), an NAD2-only matrix containing all available Parnassiana individuals and excluding all outgroups was created. Each mountain summit/population was represented by 1–7 individuals, to a total of 67 sequences. An ultrametric and bifurcating phylogenetic tree without zero branch lengths was reconstructed using BEAST 2.7 with a Birth-Death prior and 108 chain length, sampling each 1000th tree. The marker was split in two partitions to allow differences in clock rates (codons 1 + 2 and codon 3), which were calculated using a lognormal relaxed clock (Drummond et al. 2006). The site and clock models remained unlinked, while trees were linked. A gamma site model with four category counts and a GTR substitution model were applied in both partitions. Effective sample size (ESS) and stationarity were assessed for all parameters with Tracer v.1.7.1 (Rambaut et al. 2018). Results were summarized with TreeAnnotator v.2.7 (Bouckaert et al. 2019) discarding 10% of the trees as burn-in. The tree was uploaded to the GMYC web server https://species.h-its.org (Fujisawa and Barraclough 2013; Zhang 2015) and a single threshold was applied.

ABGD (Puillandre et al. 2012) and ASAP (Puillandre et al. 2021), also single-locus tools, were applied on the NAD2 and COI markers. A COI-only matrix containing all available Parnassiana individuals and excluding all outgroups was prepared. Each mountain summit/population was represented by 1–3 individuals, to a total of 47 sequences. For the NAD2-only matrix, each mountain summit/population was represented by 1–7 individuals, to a total of 67 sequences. The tests were performed in the SpartExplorer web platform (Miralles et al. 2022). Kimura (K80) was used as a substitution model and for minimum gap length different runs were tested from X = 0.5 to X = 1, as values above 1 did not provide sufficient resolution. Other parameters were set as default. Different partitions were compared and visualized through LIMES (Ducasse et al. 2020).

2.3. Morphology

A small set of morphological characters usually used as taxonomically-informative in Platycleidini was examined to complement the molecular analyses. We used these characters to assist with taxa delimitation by differentiating the clades recovered in the molecular analyses and to trace the morphological evolution in a comparative phylogenetic context. Four taxonomically significant characters used in species descriptions are compared here: for males, the last abdominal tergite (dorsal view), the cerci (dorsal view) and the titillator (anterior view); for females, the subgenital plate (ventral view). Morphological terms follow Willemse (1984) and Heller (2006). We consider the term “epiphallus” and “epiphallic sclerites”, that were used in the descriptions of species (Willemse 1973, 1980; Willemse and Willemse 1987; Heller and Willemse 1989), as synonymous to “titillator”.

The Parnassiana specimens collected during the 2021-2024 field trips were preserved in ethanol. For examination they were superficially air-dried, or, for smaller structures, submerged in ethanol, and photographed with a Zeiss Stemi 2000-C Stereo Microscope equipped with a Canon EOS 1300D camera. Genitalia were extracted and, if necessary, treated with a hot 10% KOH solution for a few minutes to remove the soft tissues. The dried specimens studied in Naturalis were placed in humid conditions overnight so the tissues would soften.

2.4. Bioacoustics

Song recordings were made in captivity in the field or in laboratory using a Tascam DR-680MKII (192 kHz, 24 bit) or, in a few occasions, ZOOM H2 (96kHz/24 bit) digital recorders, connected to a Pettersson D500 external microphone, equipped with a custom-made preamplifier. Recorded specimens were usually kept within plastic or paper containers with a volume of 400 ml covered with gauze caps. Microphones were fixed directed to the cages openings at a distance to avoid oversampling. Additional recordings were downloaded from www.xeno-canto.org. Details on specimen collecting sites, recording conditions, and equipment are provided in File S3. Recordings were visualised and analysed in Audacity v. 3.7.7 (https://www.audacityteam.org).

Bioacoustics terminology follows Heller (2006) and Ivković et al. (2017). Calling song—the song produced by an isolated male; syllable—the sound produced by one complete opening and closing movement of the tegmina; microsyllable—a short (impulse-like) syllable; macrosyllable—a long (typical) syllable; echeme—a first-order assemblage of syllables; sequence—syllable/echeme series of different length; phrase—combination of different echemes.

3. Results

3.1. Data characteristics

NAD2 was amplified for a total of 91 specimens, COI for 61 specimens, and ITS for 86 specimens. The dataset includes outgroups, all valid Parnassiana species, 27 out of 33 known Parnassiana localities reported in Willemse et al. (2018), and one new locality for the genus (Mt. Koziakas). Details and GenBank accession numbers are given in File S1. Of the three genetic markers, NAD2 has the highest variability, and the nuclear ITS the lowest. Summary statistics for all datasets are shown in Table 1, and single-locus BI trees of the three markers are shown in Files S4–S6. No signs of significant saturation or NUMTs (for the coding fragments) were found in any locus, and ESS values were above 800 for all runs. Single-locus matrices resulted in poorly-supported trees for ITS and COI, although COI mostly recovered the topology of the concatenated matrix, albeit with low support values. The NAD2 trees showed strong support with the exception of the deeper clades, and topology almost identical to the concatenated matrix. As a result of the poorly supported ITS phylogeny and the convergent topology of the strongly-supported mitochondrial markers, the trees based on the concatenated matrices reflect mostly the mitochondrial phylogeny of Parnassiana.

Table 1.
Download as
CSV
XLSX

Summary statistics for the three genetic markers and the concatenated matrix.

Genetic marker Number of individuals Sequence length (bp) Conserved sites (bp) Variable sites (bp) Parsimony informative sites (bp) Percentage of missing data
NAD2 91 993 446 547 492 2%
COI 61 969 631 338 304 34.4%
ITS 86 901 674 215 138 7.5%
Concatenated NAD2–COI–ITS matrix 93 2863 1751 1100 934

The amplified sequences were combined in a concatenated NAD2–COI–ITS matrix including all available Platycleidini outgroups and Parnassiana individuals with ≥ 2 high-quality genetic markers, which was used for ML and BI analyses and showed high ESS values (above 800 for all runs). The matrix consisted of 69 Parnassiana sequences (ingroup) and 24 outgroup sequences, and had a length of 2863 bp. For divergence time estimation, a concatenated NAD2–COI–ITS matrix which included a single individual from each mountain summit/population was preferred, leading to a total of 45 sequences (31 ingroup and 14 outgroups), and a length of 2862 bp. ESS values for all statistics were high (above 200). For species delimitation, single-locus NAD2 and COI matrices containing only Parnassiana sequences were used, with a size of 67 sequences and length of 993 bp for NAD2, and 47 sequences with a length of 969 bp for COI.

3.2. Molecular phylogeny

The phylogenetic analyses based on the concatenated matrix (Fig. 2; Files S7, S8) supported the monophyly of the genus Parnassiana. The phylogenetic tree inferred by the BI analysis (Fig. 2) showed higher support and slightly different topology than the ML tree (File S8). The ML tree provides lower support for a monophyletic Parnassiana, and thus accounts for a possible polytomy at the root of the MRCA of Parnassiana and part of Montana. As this suggests a paraphyletic Parnassiana, which contradicts the hypothesis provided by the BI phylogeny, morphology and acoustics, we prefer to base our phylogenetic inferences on the BI tree.

Based on the phylogenetic analyses, three main geographically outlined clades are revealed within Parnassiana: (1) the populations from Peloponnese, except for Mt. Taygetos (P. fusca), henceforth referred to as the Peloponnese clade; (2) the population from Mt. Taygetos, henceforth referred to as the Taygetos clade; and (3) the populations of mainland Greece, henceforth referred to as the Mainland clade (Fig. 3).

The Peloponnese clade is a highly supported monophyletic group consisting of the studied Parnassiana populations from all Peloponnesian summits except for Mt. Taygetos. The topology reflects the relative geographic position of the mountains (Fig. 3): Parnon, the southernmost and most isolated locality of the group, has a basal position; Menalon and Kyllini, in the spatial center of the group’s distribution, form a strongly supported clade with long branches; and the three mountains northwest of the former (Chelmos, Erymathos, Panachaikon) are grouped together.

The position of P. fusca, forming the Taygetos clade, is debatable, as it has low support values (compare Fig. 2 with File S8), but in no analysis it grouped with the other Peloponnesian lineages.

The well supported Mainland clade includes all lineages inhabiting Central Greece, the Pindos range, and Evvoia island (Fig. 3). It consists of four subclades whose relationships were not fully resolved due to low support at the basal nodes: (1) the Akarnanika subclade (P. nigromarginata); (2) the Dirphys-Parnassos subclade (including P. dirphys and P. parnassica); (3) the South Pindos subclade (SPS) – a well-supported lineage that spreads across most mountains of Central Greece (the southern extensions of the Pindos range); and (4) the North-Central Pindos subclade (NCP) – a group inhabiting most mountains of the Central and North Pindos range. The last two further diversify to multiple independent lineages that reveal complex geographical patterns.

The South Pindos subclade (SPS) is mainly distributed in Central Greece and consists of the currently described taxa P. gionica, P. coracis, P. tymphrestos and P. panaetolikon, as well as multiple populations with unknown taxonomic status. P. gionica, known only from Mt. Giona, branches off first. Mount Vardousia, the type locality of P. coracis, and Mt. Oxia, form the next highly supported cluster. Descending from these nodes, there is a monophyletic group including eight mountains. Mount Tymphrestos, the type locality of P. tymphrestos, has its lineage most closely related to the population of Mt. Kaliakouda. The lineage of Mt. Oiti, currently assigned to P. tymphrestos (Willemse 1973), is also related to these, though with lower support. The other branch consists of the lineage from Mt. Panaetolikon, type locality of P. panaetolikon, and a geographically distant but strongly supported cluster of lineages originating from four different mountains (Avgo, Agrafa, Gavrogo, Tzoumerka) (Fig. 3). These four lineages are scattered throughout the summits of the Central Pindos range, very distant geographically from their sister group P. panaetolikon and from all other lineages of the SPS. Despite the distant and disjunct locations of these panaetolikon-like populations (Heller and Willemse 1989; Willemse and Willemse 2008), the genetic distances between them are remarkably low. To sum up, the SPS includes almost all lineages found in the main massifs of Central Greece, as well as a taxon which has spread over a number of summits of the Central Pindos range (Parnassiana sp. 3) (Figs 2, 3).

The North-Central Pindos subclade (NCP) is the least known Parnassiana subclade, as it includes only two lineages affiliated to described species, P. tymphiensis and P. tenuis, and seven lineages that are not ascribed to taxa (compare Willemse et al. 2018). Its monophyly is well-supported, but the deeper relationships between its lineages are not resolved (Fig. 2). A well-outlined cluster within NCP is [(Tymphi + Koziakas) + (Avgo + (Karava + Voutsikaki))]. The other lineages comprising the NCP include five populations with unresolved relationships, yet clearly distinct from each other. These include: P. tenuis from Mt. Tzoumerka (type locality); a population from Mt. Voutsikaki (Parnassiana sp. 5), which is distinct from the individuals from Voutsikaki and Karava presented here as Parnassiana sp. 2; the population from Mt. Smolikas, formerly assigned to P. tymphiensis (Willemse 1980); and the lineages found on Mt. Delidimi and Mt. Helidona. Regarding the spatial distribution of NCP, all lineages but one are found in the Central and North Pindos range, with Mt. Delidimi being the southernmost, and Mt. Smolikas the northernmost locality (Fig. 3). A single lineage, however, inhabits Mt. Helidona in South Pindos, which is in close proximity to the mountains of Panaetolikon and Kaliakouda. Interestingly enough, the nuclear marker (ITS) alone places this population as sister to P. panaetolikon with high support, showing an important discrepancy between the mitochondrial and nuclear markers (see File S6).

3.3. Estimation of divergence times

The results of two of the calibration schemes – the biogeographical calibration and the NAD2 evolutionary rate for micropterous Orthopterans (Chang et al. 2020) – converged, giving estimates of less than 0.1 million years difference for the deepest nodes and much less for the more recent ones. They are presented on the BEAST chronogram in Fig. 4, while the full tree based on the micropterous molecular rate, including outgroups, can be found in the File S9. Their convergence, in addition to the strong correlation of the dating with climatic events in the area, suggest the resulting dating as the most likely scenario. Posterior probabilities and 95% HPD intervals are shown in File S10.

Figure 4. 

Time-calibrated tree of the genus Parnassiana, based on the combined NAD2–COI–ITS matrix. Numbers in white circles correspond to the estimated 95% HPD intervals and posterior probabilities in File S10. Additional node labels indicate mean node ages from two calibration schemes (above – Peloponnese separation; below – NAD2 molecular rate for micropterous Orthopterans (0.014256 subs/s/my)). The node used as a calibration point for the timescale based on the separation of the Peloponnese, dated to 4–3.5 million years ago, is marked with a black arrow. Colors of terminal nodes indicate the phylogenetic lineages and putative taxa, and correspond to the color scheme in Fig. 2. The horizontal axis represents the timescale in millions of years. Arrows below the scale refer to color-marked periods; ochre – establishment of Mediterranean climate; green – the border between Plio- and Pleistocene; light-blue – the Mid-Pleistocene Transition; dark blue – the Mid-Brunhes Transition.

The estimation based on the evolutionary rate for Pholidopterini (Çiplak et al. 2022) provided younger age estimates than the other methods of up to one million years difference for the deepest Parnassiana nodes. The resulting dates also showed low correlation to climatic events. Therefore, it is treated here as a less likely scenario and will not be further discussed. Time estimates, HPD intervals and posterior probabilities based on the Pholidopterini molecular rate can be found in File S10.

The divergence of Parnassiana from its most closely related lineage of the genus Montana is dated around 4.86/4.76 mya (evolutionary rate/biogeographical calibration) (node 1). The Most Recent Common Ancestor (MRCA) of Parnassiana, marking also the first major split into Peloponnese clade and [Taygetos clade + Mainland clade], is dated around 4.16/4.08 mya (node 2). It is followed closely by the split between the Taygetos and Mainland clades, which was used as the calibration point for the biogeographical calibration, at 3.86/3.79 mya (node 3). The poorly-supported position of P. fusca adds a measure of uncertainty at the timing of these basal nodes, but does not significantly affect the time estimates of the descending nodes, as was found by trials positioning P. fusca as a basal group of Parnassiana (trees not shown).

Descending into the main clades and subclades of Parnassiana, a number of temporal patterns can be observed. The MRCA of the Peloponnese clade (3.06/3.01 mya) and the Mainland clade (3.16/3.1 mya) were dated closely. Inside the Mainland clade, the speciation events that lead to the formation of the four subclades took place in quick succession at the end of the Pliocene, from 3.16/3.1 to 2.86/2.8 mya (nodes 4, 10, 11, 12). A second rapid radiation followed in the Early Pleistocene (2.64/2.59–2.45/2.41 mya), giving rise to the main nodes of NCP. The SPS, on the other hand, diversified later (MRCA was dated 2.24/2.19 mya; node 16), following a more gradual pattern. In total, eight out of a total of 31 cladogenic events, including major clades and subclades, took place during the Late Pliocene–Early Pleistocene (3.16–2.41 mya).

Moreover, time estimates for the majority of the cladogenic events of Parnassiana (12 nodes out of 31) fall within the Early Pleistocene (2.42/2.38 to 0.98/0.96 mya), while only a few (eight nodes, four of which correspond to a single taxon), are of Middle Pleistocene age or younger (0.77/0.75 mya). Lastly, only three lineage clusters (nodes 27–31) show time estimates matching or younger than the Mid-Brunhes Event (0.43 mya). The later data indicate long isolation periods (over one million years) between most Parnassiana populations, while only eight out of the 28 populations studied have had mitochondrial gene exchange with each other in the last 500.000 years.

3.4. Sequence-based species delimitation

Sequence-based species delimitation tests, which estimate the gap between intra- and inter-specific distances, were performed on single-locus sequences. ITS was not used due to its poor performance in resolving the phylogenetic relationships of Parnassiana.

A GMYC analysis was performed only on the NAD2 genetic marker, as it was the one with the most sequences and best bootstrap support values. The single-threshold GMYC test results showed a total of 27 ML entities with a confidence interval of 24–31 (File S12). Multiple thresholds tests overestimated the number of species significantly (results not shown).

The ABGD and ASAP methods were tested both on the NAD2 and COI markers independently. Full results are shown in the File S13. All runs of ABGD with X between 0.5–1 on NAD2 showed identical results and suggested 23 putative species, providing the least number of species across all used delimitation methods (Fig. 2). The five partitions with the best ASAP score suggest between 23 and 30 taxa. For the COI marker, all runs in ABGD with X between 0.5–1 showed identical results and suggested 18 species. Excluding the five taxa not presented in the COI dataset, the result corresponds to the ABGD results for NAD2, while for ASAP, the five partitions with the best score suggest between 20 and 25 species excluding the five taxa that were not represented.

3.5. Morphology

The comparison of morphological characters used as diagnostic in species descriptions, including male cerci and titillators (epiphallic sclerites) (Fig. 2), the shape of male 10th abdominal tergite (Fig. 5), and the female subgenital plate (Fig. 6), revealed that most terminal lineages lacked exclusive morphological synapomorphies, with only major clades recovering common patterns of their titillators and cerci. Τhe taxa belonging to the Peloponnese clade share large titillators with strong and long apical arms shaped in a similar way, and cerci that are stout and wide basally. These characters are mostly absent from the Mainland clade, with the exception of the most basal lineages representing P. nigromarginata (large titillators with strong apical arms and stout cerci) and P. dirphys and P. parnassica (stout and wide basally cerci). P. fusca exhibits intermediate characters, with the 10th abdominal tergite and cerci being similar to the Peloponnese clade, specifically to P. chelmos, while the smaller titillators lacking dark coloration with shorter apical arms seem to be transient to the Mainland clade.

Figure 5. 

Last abdominal tergite of male Parnassiana specimens from each available locality. Each locality is indicated by the mountain’s name. In case of two species co-existing in a single mountain, the taxon’s putative name is included.

Figure 6. 

Female subgenital plate of Parnassiana specimens from each available locality. Each locality is indicated by the mountain’s name. In case of two species co-existing in a single mountain, the taxon’s putative name is included.

In-group differences between the members of each major clade are less pronounced. Especially NCP shows a lot of homogeneity in the morphological characters of its members, highlighting cryptic diversity revealed by molecular data. On the other hand, there are a few striking examples of apomorphies both in basal and terminal lineages of the Mainland clade in the shape of cerci (P. nigromarginata, P. tenuis, Parnassiana sp. 2 + sp. 4), titillator (P. nigromarginata, P. tenuis, Parnassiana sp. 6), male 10th tergite and female subgenital plate (P. tenuis). Pronotum shape and relative length of tegmina also characterize some of the clades and certain lineages (compare Fig. 1). The pronotum disk is smoother and roundish in the basal lineages (Peloponnese and Taygetos clade, Akarnanika and Dirphys-Parnassos subclades). Length of tegmina relative to pronotum is shorter again in the basal lineages except for the Taygetos clade; in the SPS tegmina are longer (in males about the length of pronotum), while in the NCP male tegmina are longer than pronotum and the pronotal disk is distinctly indented from both sides of the medial keel. Dark coloration is typical for the Peloponnese clade.

3.6. Bioacoustics

In the present study we examined the song pattern in almost all known populations and phylogenetic units except for those from Oiti and Helidona. Songs are characterized by phrases consisting of series (echemes) formed by two types of syllables – microsyllable series (consisting of short syllables produced by faster movement of tegmina) and macrosyllable series (consisting of long syllables produced by slower movement of tegmina) (Figs 7, 8). Both microsyllables and macrosyllables are complex as they are produced in pairs by combined lower-amplitude faster and higher-amplitude slower movement. Song characteristics for each putative taxon are provided in File S11.

Figure 7. 

Oscillograms of male calling song examples from Parnassiana populations. Ambient temperature during recordings is noted next to the oscillograms. A1–K1, J3, K3 15-s frame; A2–K2 (except H2), D3, F3, I3, J4, J5, K4 1-s frame; H2 2-s frame; E2, F2, G2, H2, I2, I3 whole echeme; A2, B2, C2, D2, J2, J4, K2, K4 beginning of an echeme; D3, F3, J5 last part of an echeme. A P. parnon, Parnon Mt., IBER445, nocturnal recording; B P. menalon, Menalon Mt., IBER419, nocturnal recording; C P. chelmos unicolor Kyllini Mt., IBER416, diurnal recording; D P. chelmos chelmos, Chelmos Mt., IBER345, diurnal recording; E P. chelmos unicolor, Panachaikon Mt., IBER436, nocturnal recording; F P. chelmos deplanata, Erymanthos Mt., IBER369, nocturnal recording; G P. fusca, Taygetos Mt., IBER449, nocturnal recording; H P. nigromarginata, Akarnanika Mt., XC886903 (https://xeno-canto.org/886903), diurnal recording; I P. dirphys, Dirphys Mt., IBER352, diurnal recording; J1, J2 P. parnassica, Parnassos Mt., XC786736 (https://xeno-canto.org/786736), diurnal recording; J3–J5 P. parnassica, Parnassos Mt., IBER444, nocturnal recording; K1, K2 P. gionica, Giona Mt., IBER370, diurnal recording; K3, K4 P. gionica, Giona Mt., IBER370, nocturnal recording.

Figure 8. 

Oscillograms of male calling song examples from Parnassiana populations. Ambient temperature during recordings is noted next to the oscillograms. A1–N1, C2, D2, E2, J3, M2 15-s frame; A2, B2, B3, E3, E4, F2, H2–L2, J4, M3, N2 1-s frame; C3, D3, G2 2-s frame; A2, C3, D3, E3, E4, F2, G2, H2, I2, J2, J4, K2, M3, N2 whole echeme; B2 beginning of an echeme; B3 last part of an echeme and beginning of the next; L2 last part of an echeme. A P. coracis, Vardousia Mt., IBER485, nocturnal recording; B P. coracis, Vardousia Mt., IBER485, diurnal recording; C P. coracis, Oxia Mt., IBER432, diurnal recording; D1 P. tymphrestos, Tymphrestos Mt., XC886923 (https://xeno-canto.org/886923), diurnal recording; D2 P. tymphrestos, Tymphrestos Mt., XC886925 (https://xeno-canto.org/886925), nocturnal recording; E Parnassiana sp. 1, Kaliakouda Mt., IBER377, diurnal recording; F P. panaetolikon, Panaetolikon Mt., IBER440, diurnal recording; G Parnassiana sp. 3, Avgo Mt., IBER320, diurnal recording; H Parnassiana sp. 4, Avgo Mt., IBER323, diurnal recording; I Parnassiana sp. 2, Karava Mt., IBER381, diurnal recording; J1, J2 P. tymphiensis, Koziakas Mt., IBER393; J3, J4 P. tymphiensis, Koziakas Mt., IBER409; K P. tenuis, Anatolika Tzoumerka Mt., IBER313; L P. cf. tymphiensis, Smolikas Mt., IBER448, nocturnal recording; M1 Parnassiana sp. 7, Delidimi Mt., IBER349, nocturnal recording; M2, M3 Parnassiana sp. 7, Delidimi Mt., IBER349, diurnal recording; N Parnassiana sp. 5, Voutsikaki Mt., IBER490, diurnal recording.

Phylogenetic lineages within Parnassiana are characterized by the position of the microsyllable series in relation to the macrosyllable series (echemes) (before, after or isolated), and the length of the micro- and macrosyllable series defined mainly by the number of syllables within each echeme. Similarly to the morphological characterization, the basal lineages (Peloponnesian and Taygetos clade, Akarnanika, Dirphys-Parnassos subclades, but also P. gionica) group acoustically by the ancestral positioning of the microsyllable series before the macrosyllable echemes; however, in the lineages from Chelmos, Panachaiko and Erymanthos, as well as possibly in P. parnassica, the position was reversed. Uniquely, in P. nigromarginata, microsyllable series may appear both before and after the macrosyllable echemes. In the studied terminal lineages of the SPS and NCP subclades, apart from P. gionica, microsyllables are produced only after the macrosyllable echemes. In agreement with Heller’s (2006) statement, both SPS and NCP show significant variation in the duration (number of syllables) of the macrosyllable echemes, that may occur within the same individual.

4. Discussion

4.1. Phylogenetic signal in Parnassiana

The present study is a first attempt to reconstruct the phylogenetic relationships of the genus Parnassiana based on a combined set of one nuclear and two mitochondrial DNA fragments.

ITS has been previously used in Orthoptera as a nuclear marker for species and genus-level phylogenetic reconstructions (Kaya et al. 2013; Çiplak et al. 2015, 2020; Borissov et al. 2021, 2023; Kociński et al. 2022), instead of the more commonly used Histone 3, Wingless, 28S and 18S rDNA, which have much lower mutation rates and are usually preferred for higher-level phylogenies (Simon et al. 2006; Song et al. 2015, 2018; Mugleston et al. 2018). In Parnassiana, however, extensive deletions and insertions in ITS1 and ITS2, combined with a low number of parsimony informative sites, have led to low resolution and large polytomies from which few insights can be drawn, and thus the influence of the ITS matrix on the concatenated tree was limited (File S6). Similarly, low resolution of ITS for species-level phylogenies was observed in other Tettigoniinae from the tribe Pholidopterini (Çiplak et al. 2020, 2022), contrary to its high phylogenetic performance in Phaneropterinae as shown for Barbitistini (Chobanov et al. 2017; Borissov et al. 2023).

The mitochondrial marker NAD2 provided the highest resolution and statistical support among the markers (File S4), and had the most significant influence on the topology of the concatenated gene tree. NΑD2 is known to be excellent for revealing species-level phylogenies in Orthoptera, as it has a high proportion of variable and informative sites and provides well-resolved trees (Simon et al. 1994, 2006; Çiplak et al. 2015, 2020; Borissov et al. 2021; Kociński et al. 2022; Willemse et al. 2023). COI, a marker very common in comparative and species-delimitation studies, provided similar topology to NAD2, albeit with fewer parsimony-informative and variable sites and lower support values (File S5). Such results indicating that NAD2 has more parsimony-informative and variable sites than COI and higher branch support are in line with previous observations in Orthoptera genera, such as Poecilimon and Isophya (Chobanov et al. 2017; Borissov and Chobanov 2020), and in Odonata (Cheng et al. 2018).

Comparison of the mitochondrial and nuclear phylogenies can reveal signs of mitonuclear discordance. Due to the poor resolution of the ITS tree, however, the concatenated gene tree largely reflects the mitochondrial evolutionary history of Parnassiana, and we rely on signs from phenological characters, mostly those subjected to sexual selection, i.e., genitalia and acoustic communication, to detect such discordances.

4.2. Internal phylogeny of Parnassiana and its reflection on the evolution of morphological and acoustic traits

The phylogenetic tree presented here is mostly congruent with the current systematics of the genus (Massa and Fontana 2011), as it supports the monophyly of Parnassiana and all thirteen nominal species. Furthermore, the phylogenetic analysis and species delimitation tests support the observation by Willemse and Willemse (2008) that certain populations found in isolated massifs constitute unnamed species.

Signs of discrepancies between the presented gene tree and species delineations based on morphology are observed in several cases. Τhere is at least one clear example of mitochondrial introgression in our phylogeny. Within the SPS, Parnassiana sp. 1 groups with P. tymphrestos but shows clear morphological and acoustic characteristics of the P. panaetolikon + Parnassiana sp. 3 group. In the Peloponnese clade, the lineage from Mt. Kyllini (P. chelmos unicolor), formerly grouped under the same taxon as the Panachaikon population (Willemse 1973, 1980), is sister to P. menalon instead of the other populations currently placed within P. chelmos (Willemse 1973, 1975, 1980). The molecular phylogeny is further supported by the song pattern of the Kyllini population (Fig. 7C), which clearly groups with P. parnon and P. menalon based on the long phrases starting with microsyllable series, while in other members of P. chelmos the microsyllable series are produced at the end of the phrases (Fig. 7D–F). Yet, the Kyllini lineage shows clear morphological similarity with P. chelmos (Figs 1, 2, 5, 6). There are two possible scenarios to explain this pattern. One scenario is that the Kyllini population kept mitochondrial DNA from an extinct sister species to P. menalon that once existed on Kyllini, whose gene pool was flooded by immigrating males from the surrounding mountains hosting P. chelmos, a process called ghost introgression (Shen et al. 2025). This scenario requires a convergent evolution of the song pattern in Kyllini males that may be driven by sexual selection of the females of the local species. Another explanation of the observed phenomenon may be that the morphological characters of P. chelmos are plesiomorphies shared by the entire clade, whereas the unique morphology of P. menalon presents an apomorphy or atavistic reappearance (a scenario that is partly supported by the ITS phylogeny – File S6). Another intriguing case which could be related to ghost introgression is Parnassiana sp. 6 (Mt. Helidona). Based on the mitochondrial phylogeny, the population belongs to the NCP, but geographically it belongs in Central Greece and shows morphological similarity with the P. panaetolikon and Parnassiana sp. 1 from its neighboring mountains Panaetolikon and Kaliakouda. Moreover, the ITS data clearly group it with P. panaetolikon.

Such discrepancies suggest that the mechanisms that determine the evolution of Parnassiana are complex, and a simple allopatric speciation model may not be sufficient to interpret them, a pattern that is common among sky island inhabitants (Recuero et al. 2014; Knowles and Massatti 2017; Ortego and Knowles 2022).

Based on such an interesting phylogenetic pattern, formerly discussed phenological peculiarities in the genus may be discussed in an evolutionary context. Former studies supported the view that speciation processes in Parnassiana are significantly influenced by sexual selection on male genitalia (titillators, last tergite, cerci) and less on the evolution of calling songs, which are similar between most described taxa and variable within species, therefore having poor systematic value (Heller 2006). Indeed, morphology characterizes the main phylogenetic lineages, with the Peloponnese and Taygetos clades, and the Akarnanika and Dirphys-Parnassos subclades showing clear morphological differentiation. Additionally, a few outliers of the mostly uniform morphology with two types of titillators (short- and long-armed) observed at the terminal nodes of the SPS and NCP, can be observed in P. tenuis (various characters) and Parnassiana sp. 4 (in cerci), both sympatric with Parnassiana sp. 3.

Song characteristics also support the main lineages, though their pattern is different from the one observed in the morphology, suggesting a distinct evolution from it. Here, main song-type groups outline the following groupings: 1. Parnon+Menalo+Kyllini+Taygetos (long to medium macrosyllable echemes with microsyllable series produced isolated or before the macrosyllable series); 2. Chelmos+Panachaiko+Erymanthos (short to medium-length macrosyllable echemes with long microsyllable series produced after the macrosyllable series); 3. Akarnanika (unique pattern with microsyllable series produced both before and after the macrosyllable series); 4. Dirphys (short macrosyllable echemes with microsyllable series produced before the macrosyllable series containing long macrosyllables); 5. Parnassos and Giona (two types of songs with very long nocturnal macrosyllable echemes but differing in the production of microsyllables); and 6. the rest of the SPS together with NCP (variable short to medium-length macrosyllable echemes with short microsyllable series produced after the macrosyllable series). It is worth noting the resemblance of song types of distantly related lineages like P. chelmos, the SPS (except P. gionica), and the NCP, that obviously evolved independently superficially similar song pattern apomorphies.

Significant intraindividual variation was observed within the SPS and NCP subclades depending on the ambient temperature, with specimens tending to produce denser and shorter (composed of less syllables) macrosyllable echemes with increased temperature (examples from own recordings with known details shown in Fig. 8J1 and J3, M1 and M2). At the same time, continuously singing animals tend to change the temporal song pattern with time, which may reflect change in the body temperature (warming up) at similar ambient temperature (Fig. 8C1 and C2, E1 and E2). Yet, these variations do not reach the scope of differentiation in song patterns between sympatric taxa, i.e., Parnassiana sp. 3 (Fig. 8G) occurring together with Parnassiana sp. 4 (Fig. 8H) and P. tenuis (Fig. 8K).

Regarding the sympatric taxa, Heller and Willemse (1989) suggested that the songs of two taxa found in the same location (P. tenuis and an unnamed species in Tzoumerka) show very limited differences in song pattern compared to their morphological distinction. However, according to our observations, the P. tenuis song considerably differs from the syntopic Parnassiana sp. 3, having much shorter macrosyllable echemes (compare Figs 8G and K). Similar song differentiation was observed in the syntopic Parnassiana sp. 3 and sp. 4 on Mt. Avgo (Figs 8G and H).

Much more interesting are the observed examples of different echeme length depending on the time of the day. Although there are no distinct diurnal and nocturnal song patterns in Parnassiana (contrary to some related groups like Montana; see Ivković et al. 2017), the number of syllables within the macrosyllable echemes and thus their length may differ when songs are produced in day or night (see Fig. 7J, K for examples of diurnal and nocturnal songs produced at similar temperature). In P. parnassica and P. gionica, nocturnal songs (Fig. 7J3, K3) contain many more macrosyllables within the series than diurnal songs (Fig. 7J1, K1) at similar temperature. In those cases, temperature seems not to be a major factor and thus the case may concern a primitive stage of distinction between diurnal and nocturnal song patterns known in some other related genera (e.g. Ivković et al. 2017; Barataud 2025).

The above discussed phenotypic evolution of morphological and acoustic traits of Parnassiana supports the complex phylogenetic pattern and proves that speciation processes in the group have been led by sexual selection that accelerated diversification in sympatric taxa. Both the evolution of genitalia and behavior involved in reproduction show strong phylogenetic signals but followed distinct evolutionary pathways.

4.3. Evolution of Parnassiana in a paleogeographic and paleoclimatic context

Parnassiana originated in the southern areas of the Balkan Peninsula ca. 4.86/4.76 mya in the Middle Pliocene from a common ancestor shared with a yet poorly outlined infragroup of Montana (Fig. 4). Shortly after, the Corinthian Rift developed, causing the isolation of Peloponnese from mainland Greece 4–3.5 mya (Wicker et al. 2024), and the first major split of Parnassiana into Peloponnesian clade + [Taygetos + Mainland clades] (4.16–4.08 mya). Around the same time (4–2 mya) the uplifting of Mt. Taygetos occurred, which, together with rifting that formed deep and occasionally flooded valleys, led to the massif’s isolation from the rest of the Peloponnese (Fountoulis 2014; Kleman et al. 2016) and to the permanent isolation of P. fusca both from the highly-supported Peloponnesian clade and from the Mainland clade. These two processes resulted in the geographic split of the ancestral lineages of Parnassiana and the establishment of its three main clades (Peloponnesian, Taygetos and Mainland) via vicariance. This timeline of divergence 4 to 3 mya is shared by many different organisms in the area, such as the land snail Codringtonia (Kotsakiozi et al. 2012), the mammal Talpa stankovici (Colangelo et al. 2010), lizards (Psonis et al. 2018), the scorpion genus Euscorpius (Parmakelis et al. 2013), and a number of Dolichopoda species (Allegrucci et al. 2009, 2021).

In the Late Pliocene (3.6–2.58 mya; Jiménez-Moreno et al. 2013), constant climate cooling and drying transformed the savannas and subtropical humid forests of the area to dominating sclerophyllous and steppe vegetation with the establishment of drier seasonal Mediterranean climate 3.2 to 2.8 mya (Suc 1984). This process has been shown to cause extensive diversification across the Mediterranean in Orthoptera (Çiplak et al. 2015; Allegrucci et al. 2021; Uluar et al. 2023) and other organisms (Kotsakiozi et al. 2012; Fiz-Palacios and Valcárcel 2013; Poulakakis et al. 2015; Kougioumoutzis et al. 2021; Koutroumpa et al. 2021). This cooler and drier climate would have promoted wider establishment of the cold-adapted Parnassiana in the high altitude mountain ranges of the southern Balkan Peninsula. At the same time, long-distance dispersal and contact between populations was prevented due to the fragmentation of the most isolated mountain massifs by a network of lowlands with unsuitable climate. As a result, during this period, migration and isolation events may have promoted secondary diversification within the Peloponnese clade and the establishment of the four major subclades within the Mainland clade. In fact, eight out of a total of 31 cladogenic events, including major clades and subclades, took place during the Late Pliocene–Early Pleistocene (3.16–2.41 mya), marking this time period as a major time stamp in the evolution of Parnassiana, during which the major clades and subclades of the genus were established and the NCP showed a burst of early diversification. In the Pleistocene, climate continued to change and climate cycles that marked the evolutionary pattern of the European biota were established (Hewitt 1996, 2004; Taberlet 1998; Wallis et al. 2016). During the Early Pleistocene (2.58–0.77 mya), short glacial periods of low intensity alternated with interglacials at 41 kyr cycles (Willeit et al. 2019; Shackleton et al. 2023). The majority of extant Parnassiana lineages evolved during this period (Fig. 4), possibly as a result of short stepping stone dispersals during the mild glacial periods followed by isolation periods preventing gene flow between the populations due to the warm interglacials. By the Mid-Pleistocene Transition 1.25–0.7 mya, which marks the end of this stage (Willeit et al. 2019; Shackleton et al. 2023), all putative Parnassiana species were established.

The Mid-Pleistocene Transition marked a significant change in Earth’s climate from “mild” to “full” glacial periods (Willeit et al. 2019; Shackleton et al. 2023). Glacial-interglacial cycles became longer, from 41.000-year to 100.000-year periods, with shorter interglacials, and longer and more intense glacial periods (Barth et al. 2018; Shackleton et al. 2023). Nevertheless, the intense glacial-interglacial cycles of the Middle and Late Pleistocene appear to have had limited effects on the speciation of Parnassiana. Only a small number of lineages diversified after the Mid-Pleistocene Transition (e.g., P. coracis, P. chelmos), and these splits produce only intraspecific divergence. Only three taxa show signs of recent long-distance dispersal (either genetic and/or geographic – Parnassiana sp. 3, P. tymphiensis and P. tenuis), taking advantage of the intense glacial periods to expand their ranges using the Pindos mountains as a corridor.

The evolutionary history of Parnassiana is reflected in its geographic patterns. Massifs isolated by larger distance and low-altitude valleys, such as Taygetos, Dirphys, and Akarnanika, host well-defined species with long evolutionary histories, reflected in their morphological and song apomorphies. On the other hand, summits that are parts of large massifs, such as the mountains of northern Peloponnese, the Central Greece mountains and the Pindos range, are characterized by taxa with more uniform phenotype and complex phylogenetic relationships, with some exceptions that may represent splits of currently isolated or extinct lineages. Such patterns have been previously observed for inhabitants of sky islands such as grasshoppers, spiders and birds (Masta 2000; Robin et al. 2010; He and Jiang 2014; Ortego and Knowles 2022). Secondary dispersals of Parnassiana are characteristic for the Pindos range and are facilitated by the extensive high-altitude range and thus good connectivity and cooler climate. Populations here tend to have wider distributions and higher densities than populations inhabiting isolated small mountain summits (like P. parnon, P. menalon, P. nigromarginata, P. dirphys; own observation). Such dispersals and secondary contacts resulted in the complex phylogenetic pattern of the crown clades: populations from distant mountains grouping in the tree (P. panaetolikon and Parnassiana sp. 3); neighboring mountains hosting distant lineages (Helidona and Kaliakouda); and distinct lineages with syntopic occurrence (Tzoumerka, Avgo, and Voutsikaki).

4.4. Niche conservatism and refugia in Parnassiana

The clear correlations between cladogenic events, geographical distribution and geo-climatic events, suggest that allopatric speciation was the predominant, though not exclusive, mechanism shaping the evolution and speciation of Parnassiana. This is a very common pattern observed in sky islands (He and Jiang 2014; Martinez-Sañudo et al. 2022), and is primarily driven by niche conservatism, “the tendency of lineages to maintain their ancestral ecological niche” (Wiens 2004). Niche conservatism seems to be especially important in Parnassiana (Stefanidis et al. 2024, 2025), as is the case for many montane species (He and Jiang 2014). Close relatives of Parnassiana, such as Montana and Modestana, are typical of cool climates with most species occurring in mountains or cool deserts (Massa and Fontana 2011). It is therefore expected that the ancient taxon that gave rise to Parnassiana would also share these climatic preferences, and inherit them to its descendants.

Parnassiana populations diversified during the Late Pliocene and Early Pleistocene. The old divergence times indicate that populations, with a few exceptions (Parnassiana sp. 3, P. chelmos, P tymphiensis, P. tenuis), survived locally in the sky islands that acted as interglacial refugia (see Berger et al. 2010 and Çiplak et al. 2015 for similar scenarios), without successfully colonizing nearby mountains, through several glacial cycles. A remarkably similar evolutionary pattern and history was observed in Anterastes, a bush-cricket genus with sky-island distribution in Anatolia (Çiplak et al. 2015; Uluar et al. 2023), highlighting the importance of Pliocene–Early Pleistocene events in the Orthoptera inhabiting the sky islands of the Eastern Mediterranean.

4.5. Species delimitation

Species delimitation analyses based on the NAD2 and COI molecular markers suggest a range of 18–41 taxa in Parnassiana (Fig. 2; Files S12, S13). ABGD based on the complete sample of NAD2 provided the most conservative number of species with 23 putative taxa, which is congruent with the lowest of the five best estimates of the NAD2-ASAP tests. The NAD2-based GMYC analysis, however, suggested a larger number of taxa (27). GMYC is known for overestimating the number of species, especially when dealing with a small number of individuals per species and strong within-species divergence, which is frequent in groups with geographically isolated populations like Parnassiana (Luo et al. 2018; Hofmann et al. 2019). Additionally, species delimitation analyses based on single-locus DNA do not account for gene tree vs species tree discordance, incomplete lineage sorting, introgression, or gene flow (Luo et al. 2018; Puillandre et al. 2021), processes which appear to have affected the evolution of Parnassiana. Taking into account these inherent biases and the studied phenotype characteristics, we expect that the true number of taxa falls within the lowest estimates of the delimitation tests.

A few notable discrepancies between species delimitation tests and phenotypic characters are discussed below. The ABGD test on the NAD2 marker suggested the lumping of the population of Mt. Kaliakouda with P. tymphrestos from Mt. Tymphrestos, even though the titillators of the two populations bear remarkable differences, and GMYC and ASAP results indicate that the two are distinct species (Fig. 2). As a result, we suggest that Mt. Kaliakouda hosts a yet undescribed species of Parnassiana. On the opposite end, even though all methods recover the specimens from the mountains hosting P. chelmos (excluding Mt. Kyllini) as two or even three distinct entities (Fig. 2), they can be regarded as subspecies based on their conserved phenotype, as has already been proposed (Willemse 1973; Willemse 1980). Similarly, despite GMYC and ASAP suggesting that the lineage pairs inhabiting the neighboring mountains of Vardousia and Oxia, as well as Karava and Voutsikaki, are distinct entities, we consider them as infraspecific based on the ABGD results, the recent divergence times (younger than 0.6 mya), and the lack of discernible phenotypic differentiation. Therefore, we accept the placement of the Parnassiana population of Mt. Oxia in P. coracis, as suggested by Stefanidis et al. (2025), and suggest that the lineages of Karava and Voutsikaki constitute a single undescribed species.

Overall, based on the present results including phylogenetic relationships, species delimitation tests, and a combination of phenotypic characters, we suggest the occurrence of 23 putative species of Parnassiana in Greece.

5. Conclusions

In this study, we found clear correlations between cladogenic events in Parnassiana, its geographical distribution and geo-climatic events, suggesting that allopatric speciation was the predominant, though not exclusive, mechanism shaping the evolutionary history of the genus. Significant paleogeographic and paleoclimatic events that caused major lineage splits in Parnassiana include: (1) the isolation of Peloponnese from Central Greece 4–3.5 mya; (2) the establishment of the Mediterranean climate 3.2–2.8 mya; (3) the Plio-Pleistocene Transition 2.58 mya; and (4) the glacial-interglacial cycles of Early Pliocene. All putative Parnassiana species had been established before the onset of the intense glacial cycles of the Middle and Late Pleistocene, which were instead characterized mainly by dispersal, intra-species diversification, and possible gene exchange.

Speciation in Parnassiana was ruled by sexual selection in allopatry, with morphological and acoustic traits connected with reproduction changing gradually over its evolutionary history. Major lineages maintained basic morphological characteristics and song patterns unless in secondary contact, during which distant lineages experienced reinforcement, while in other cases populations intermixed and exchanged genetic information leading to the establishment of ghost mitogenomes. A combination of species delimitation tests, morphological characters, and bioacoustics, suggest the presence of 23 putative Parnassiana species in Greece.

6. Declarations

Conflict of interest. The authors have declared that no competing interests exist.

Author contributions. All authors listed have made a substantial, direct and intellectual contribution to the work, and approved it for publication.

7. Acknowledgements

This study is part of grant KP-06-N81/5–04.12.2024 to Dragan Chobanov by the National Science Fund (MES) of Bulgaria. Additional grants that supported this study include the Theodore J. Cohn Research Fund from the Orthopterists’ Society and a Grant supporting the Orthoptera Species File (Orthoptera of the Balkan Peninsula and the Carpathian Basin II: a database of digital data in the Orthoptera Species File). The Fonts Pontium fund from Naturalis Biodiversity Center, Leiden, provided the opportunity to examine important material stored in its collections. The samples underlying this study are preserved in the facilities upgraded by project DiSSCo-BG funded by the National Roadmap for Research Infrastructures, Ministry of Education and Science of the Republic of Bulgaria.

Special thanks go to Luc Willemse, Charlotte Hartong, Anna Ruijbroek, and the Naturalis Biodiversity Center for their warm welcome and support during the visit.

Material was collected with permissions 62419/1948 from 26/7/2021 and 71234/2214 from 27/7/2023 issued by the Greek authorities (Greece). We warmly thank Apostolis Stefanidis, Vassiliki Kati, Manolis Avramakis and Luc Willemse for collecting specimens from a few localities in Central Greece (Oiti, Parnassos, Helidona).

We sincerely thank the two reviewers, Martin Husemann and Mattia Ragazzini, and the editor, Lara-Sophie Dey, for the constructive comments and suggestions, which improved this manuscript.

8. References

  • Allegrucci G, Rampini M, Gratton P, Todisco V, Sbordoni V (2009) Testing phylogenetic hypotheses for reconstructing the evolutionary history of Dolichopoda cave crickets in the eastern Mediterranean. Journal of Biogeography 36(9): 1785–1797. https://doi.org/10.1111/j.1365-2699.2009.02130.x
  • Allegrucci G, Rampini M, Chimenti C, Alexiou S, Di Russo C (2021) Dolichopoda cave crickets from Peloponnese (Orthoptera, Rhaphidophoridae): molecular and morphological investigations reveal four new species for Greece. The European Zoological Journal 88(1): 505–524. https://doi.org/10.1080/24750263.2021.1902005
  • Berger D, Chobanov DP, Mayer F (2010) Interglacial refugia and range shifts of the alpine grasshopper Stenobothrus cotticus (Orthoptera: Acrididae: Gomphocerinae). Organisms Diversity & Evolution 10(2): 123–133. https://doi.org/10.1007/s13127-010-0004-4
  • Borissov SB, Hristov GH, Chobanov DP (2021) Phylogeography of the Poecilimon ampliatus species group (Orthoptera: Tettigoniidae) in the context of the Pleistocene glacial cycles and the origin of the only thelytokous parthenogenetic phaneropterine bush-cricket. Arthropod Systematics & Phylogeny 79: 401–418. https://doi.org/10.3897/asp.79.e66319
  • Borissov SB, Bobeva A, Çıplak B, Chobanov DP (2020) Evolution of Poecilimon jonicus group (Orthoptera: Tettigoniidae): a history linked to the Aegean Neogene paleogeography. Organisms Diversity & Evolution 20(4): 803–819. https://doi.org/10.1007/s13127-020-00466-9
  • Borissov SB, Heller K, Çıplak B, Chobanov DP (2023) Origin, evolution and systematics of the genus Poecilimon (Orthoptera: Tettigoniidae)—An outburst of diversification in the Aegean area. Systematic Entomology 48(1): 198–220. https://doi.org/10.1111/syen.12580
  • Bouckaert R, Vaughan TG, Barido-Sottani J, Duchêne S, Fourment M, Gavryushkina A, Heled J, Jones G, Kühnert D, Maio ND, Matschiner M, Mendes FK, Müller NF, Ogilvie HA, Plessis L du, Popinga A, Rambaut A, Rasmussen D, Siveroni I, Suchard MA, Wu C-H, Xie D, Zhang C, Stadler T, Drummond AJ (2019) BEAST 2.5: An advanced software platform for Bayesian evolutionary analysis (v.2.7). PLOS Computational Biology 15(4): e1006650. https://doi.org/10.1371/journal.pcbi.1006650
  • Chang H, Qiu Z, Yuan H, Wang X, Li X, Sun H, Guo X, Lu Y, Feng X, Majid M, Huang Y (2020) Evolutionary rates of and selective constraints on the mitochondrial genomes of Orthoptera insects with different wing types. Molecular Phylogenetics and Evolution 145: 106734. https://doi.org/10.1016/j.ympev.2020.106734
  • Cheng Y-C, Chen M-Y, Wang J-F, Liang A-P, Lin C-P (2018) Some mitochondrial genes perform better for damselfly phylogenetics: species-and population-level analyses of four complete mitogenomes of Euphaea sibling species. Systematic Entomology 43(4): 702–715. https://doi.org/10.1111/syen.12299
  • Chernomor O, Von Haeseler A, Minh BQ (2016) Terrace aware data structure for phylogenomic inference from supermatrices. Systematic biology 65(6): 997–1008. https://doi.org/10.1093/sysbio/syw037
  • Chobanov DP, Kaya S, Grzywacz B, Warchalowska-Śliwa E, Çıplak B (2017) The Anatolio-Balkan phylogeographic fault: A snapshot from the genus Isophya (Orthoptera, Tettigoniidae). Zoologica Scripta 46(2): 165–179. https://doi.org/10.1111/zsc.12194
  • Çiplak B, Yahyaoğlu Ö, Uluar O (2020) Revisiting Pholidopterini (Orthoptera, Tettigoniidae): Rapid radiation causes homoplasy and phylogenetic instability. Zoologica Scripta 50(2): 225–240. https://doi.org/10.1111/zsc.12463
  • Çiplak B, Kaya S, Boztepe Z, Gündüz İ (2015) Mountainous genus Anterastes (Orthoptera, Tettigoniidae): autochthonous survival across several glacial ages via vertical range shifts. Zoologica Scripta 44(5): 534–549. https://doi.org/10.1111/zsc.12118
  • Çiplak B, Yahyaoğlu Ö, Uluar O, Doğan Ö, Başibüyük HH, Korkmaz EM (2022) Phylogeography of Pholidopterini: Revising molecular clock calibration by Mid-Aegean Trench. Insect Systematics & Evolution 53(5): 515–535. https://doi.org/10.1163/1876312X-bja10033
  • Colangelo P, Bannikova AA, Kryštufek B, Lebedev VS, Annesi F, Capanna E, Loy A (2010) Molecular systematics and evolutionary biogeography of the genus Talpa (Soricomorpha: Talpidae). Molecular Phylogenetics and Evolution 55(2): 372–380. https://doi.org/10.1016/j.ympev.2010.01.038
  • De Baets K, Antonelli A, Donoghue PCJ (2016) Tectonic blocks and molecular clocks. Philosophical Transactions of the Royal Society B: Biological Sciences 371(1699): 20160098. https://doi.org/10.1098/rstb.2016.0098
  • Fassoulas CG (2018) The geodynamic and paleogeographic evolution of the Aegean in the Tertiary and Quaternary: A review. In: Sfenthourakis et al. S (Ed.), Biogeography and Biodiversity of the Aegean. In honour of Prof. Moysis Mylonas. Broken Hill Publishers Ltd, Nicosia, Cyprus, 25–45.
  • Fiz-Palacios O, Valcárcel V (2013) From Messinian crisis to Mediterranean climate: A temporal gap of diversification recovered from multiple plant phylogenies. Perspectives in Plant Ecology, Evolution and Systematics 15(2): 130–137. https://doi.org/10.1016/j.ppees.2013.02.002
  • Fountoulis I (2014) Quaternary Basin Sedimentation and Geodynamics in SW Peloponnese (Greece) and Late Stage Uplift of Taygetos Mt. Bollettino di Geofisica Teorica ed Applicata 55(2): 303–324. https://doi.org/10.4430/bgta0074
  • Fujisawa T, Barraclough TG (2013) Delimiting Species Using Single-Locus Data and the Generalized Mixed Yule Coalescent Approach: A Revised Method and Evaluation on Simulated Data Sets. Systematic Biology 62(5): 707–724. https://doi.org/10.1093/sysbio/syt033
  • Gawthorpe RL, Leeder MR, Kranis H, Skourtsos E, Andrews JE, Henstra GA, Mack GH, Muravchik M, Turner JA, Stamatakis M (2018) Tectono-sedimentary evolution of the Plio-Pleistocene Corinth rift, Greece. Basin Research 30(3): 448–479. https://doi.org/10.1111/bre.12260
  • Gkontas I, Papadaki S, Trichas A, Poulakakis N (2016) First assessment on the molecular phylogeny and phylogeography of the species Gnaptor boryi distributed in Greece (Coleoptera: Tenebrionidae). Mitochondrial DNA Part A 28(6): 927–934. https://doi.org/10.1080/24701394.2016.1209196
  • Heller K-G (2006) Song evolution and speciation in bushcrickets. In: Drosopoulos S, Claridge M (Eds), Insect Sounds and Communication: Physiology, Behaviour, Ecology, and Evolution. CRC Press at Taylor & Francis Group, Boca Raton, 137–152. Available from: https://doi.org/10.1603/008.102.0420.
  • Heller K-G, Willemse F (1989) Two new bush-crickets from Greece, Leptophyes lisae sp. nov. and Platycleis (Parnassiana) tenuis sp. nov. (Orthoptera: Tettigoniidae). Entomologische Berichten 49(10): 144–156. https://natuurtijdschriften.nl/pub/1012545
  • Hemelsdaël R, Charreau J, Ford M, Proborukmi MS, Malartre F, Urban B, Blard P-H (2021) Tectono-climatic controls of the early rift alluvial succession: Plio-Pleistocene Corinth Rift (Greece). Palaeogeography, Palaeoclimatology, Palaeoecology 576: 110507. https://doi.org/10.1016/j.palaeo.2021.110507
  • Hewitt G (1996) Some genetic consequences of ice ages, and their role in divergence and speciation. Biological Journal of the Linnean Society 58(3): 247–276. https://doi.org/10.1006/bijl.1996.0035
  • Hoang DT, Chernomor O, Von Haeseler A, Minh BQ, Vinh LS (2018) UFBoot2: improving the ultrafast bootstrap approximation. Molecular biology and evolution 35(2): 518–522. https://doi.org/10.1093/molbev/msx281
  • Hofmann EP, Nicholson KE, Luque-Montes IR, Köhler G, Cerrato-Mendoza CA, Medina-Flores M, Wilson LD, Townsend JH (2019) Cryptic diversity, but to what extent? Discordance between single-locus species delimitation methods within mainland anoles (Squamata: Dactyloidae) of northern central America. Frontiers in Genetics 10: 11. https://doi.org/10.3389/fgene.2019.00011
  • Ivković S, Iorgu IS, Horvat L, Chobanov D, Korsunovskaya O, Heller K-G (2017) New data on the bush-cricket Montana medvedevi (Orthoptera: Tettigoniidae), critically endangered in Europe (EU 28), and a comparison of its song with all known song patterns within the genus. Zootaxa 4263(3): 527–542. https://doi.org/10.11646/zootaxa.4263.3.5
  • Jiménez-Moreno G, Burjachs F, Expósito I, Oms O, Carrancho Á, Villalaín JJ, Agustí J, Campeny G, Gómez De Soler B, Van Der Made J (2013) Late Pliocene vegetation and orbital-scale climate changes from the western Mediterranean area. Global and Planetary Change 108: 15–28. https://doi.org/10.1016/j.gloplacha.2013.05.012
  • Kalyaanamoorthy S, Minh BQ, Wong TK, Von Haeseler A, Jermiin LS (2017) ModelFinder: fast model selection for accurate phylogenetic estimates. Nature methods 14(6): 587–589. https://doi.org/10.1038/nmeth.4285
  • Kaya S, Boztepe Z, Çiplak B (2013) Phylogeography of Troglophilus (Orthoptera: Troglophilinae) based on Anatolian members of the genus: radiation of an old lineage following the Messinian. Biological Journal of the Linnean Society 108(2): 335–348. https://doi.org/10.1111/j.1095-8312.2012.02025.x
  • Kleman J, Borgström I, Skelton A, Hall A (2016) Landscape evolution and landform inheritance in tectonically active regions: The case of the Southwestern Peloponnese, Greece. Zeitschrift für Geomorphologie 60(2): 171–193. https://doi.org/10.1127/zfg/2016/0283
  • Knowles LL, Massatti R (2017) Distributional shifts – not geographic isolation – as a probable driver of montane species divergence. Ecography 40(12): 1475–1485. https://doi.org/10.1111/ecog.02893
  • Kociński M, Chobanov D, Grzywacz B (2022) New insights into the genetic diversity of the Balkan bush-crickets of the Poecilimon ornatus group (Orthoptera: Tettigoniidae). Arthropod Systematics & Phylogeny 80: 243–259. https://doi.org/10.3897/asp.80.e82447
  • Kornilios P, Thanou E, Kapli P, Parmakelis A, Chatzaki M (2016) Peeking through the trapdoor: Historical biogeography of the Aegean endemic spider Cyrtocarenum Ausserer, 1871 with an estimation of mtDNA substitution rates for Mygalomorphae. Molecular Phylogenetics and Evolution 98: 300–313. https://doi.org/10.1016/j.ympev.2016.01.021
  • Kotsakiozi P, Parmakelis A, Giokas S, Papanikolaou I, Valakos ED (2012) Mitochondrial phylogeny and biogeographic history of the Greek endemic land-snail genus Codringtonia Kobelt 1898 (Gastropoda, Pulmonata, Helicidae). Molecular Phylogenetics and Evolution 62(2): 681–692. https://doi.org/10.1016/j.ympev.2011.11.012
  • Kougioumoutzis K, Kokkoris I, Panitsa M, Kallimanis A, Strid A, Dimopoulos P (2021) Plant Endemism Centres and Biodiversity Hotspots in Greece. Biology 10(2): 72. https://doi.org/10.3390/biology10020072
  • Koutroumpa K, Warren BH, Theodoridis S, Coiro M, Romeiras MM, Jiménez A, Conti E (2021) Geo-Climatic Changes and Apomixis as Major Drivers of Diversification in the Mediterranean Sea Lavenders (Limonium Mill.). Frontiers in Plant Science 11: 612258. https://doi.org/10.3389/fpls.2020.612258
  • Kumar S, Stecher G, Li M, Knyaz C, Tamura K (2018) MEGA X: molecular evolutionary genetics analysis across computing platforms. Molecular Biology and Evolution 35(6): 1547–1549. https://doi.org/10.1093/molbev/msy096
  • Luo A, Ling C, Ho SY, Zhu C-D (2018) Comparison of methods for molecular species delimitation across a range of speciation scenarios. Systematic Biology 67(5): 830–846. https://doi.org/10.1093/sysbio/syy011
  • Martinez-Sañudo I, Basso A, Ortis G, Marangoni F, Stancher G, Mazzon L (2022) Strong genetic differentiation between fragmented alpine bush-cricket populations demands preservation of evolutionary significant units. Insect Conservation and Diversity 15(6): 752–762. https://doi.org/10.1111/icad.12601
  • Massa B, Fontana P (2011) Supraspecific taxonomy of Palaearctic Platycleidini with unarmed prosternum: a morphological approach (Orthoptera: Tettigoniidae, Tettigoniinae). Zootaxa 2837(1): 1–47. https://doi.org/10.11646/zootaxa.2837.1.1
  • Miralles A, Ducasse J, Brouillet S, Flouri T, Fujisawa T, Kapli P, Knowles L.L., Kumari S, Stamatakis A, Sukumaran J, Lutteropp S, Vences M. Puillandre N (2022) SPART: A versatile and standardized data exchange format for species partition information. Molecular Ecology Resources 22(1): 430–438. https://doi.org/10.1111/1755-0998.13470
  • Mugleston JD, Naegle M, Song H, Whiting MF (2018) A Comprehensive Phylogeny of Tettigoniidae (Orthoptera: Ensifera) Reveals Extensive Ecomorph Convergence and Widespread Taxonomic Incongruence. Insect Systematics and Diversity 2(4). https://doi.org/10.1093/isd/ixy010
  • Nguyen L-T, Schmidt HA, Von Haeseler A, Minh BQ (2015) IQ-TREE: a fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Molecular biology and evolution 32(1): 268–274. https://doi.org/10.1093/molbev/msu300
  • Ortego J, Knowles LL (2022) Geographical isolation versus dispersal: Relictual alpine grasshoppers support a model of interglacial diversification with limited hybridization. Molecular Ecology 31(1): 296–312. https://doi.org/10.1111/mec.16225
  • Parmakelis A, Kotsakiozi P, Stathi I, Poulikarakou S, Fet V (2013) Hidden diversity of Euscorpius (Scorpiones: Euscorpiidae) in Greece revealed by multilocus species-delimitation approaches. Biological Journal of the Linnean Society 110(4): 728–748. https://doi.org/10.1111/bij.12170
  • Poulakakis N, Kapli P, Lymberakis P, Trichas A, Vardinoyiannis K, Sfenthourakis S, Mylonas M (2015) A review of phylogeographic analyses of animal taxa from the Aegean and surrounding regions. Journal of Zoological Systematics and Evolutionary Research 53(1): 18–32. https://doi.org/10.1111/jzs.12071
  • Psonis N, Antoniou A, Karameta E, Leaché AD, Kotsakiozi P, Darriba D, Kozlov A, Stamatakis A, Poursanidis D, Kukushkin O, Jablonski D, Crnobrnja–Isailović J, Gherghel I, Lymberakis P, Poulakakis N (2018) Resolving complex phylogeographic patterns in the Balkan Peninsula using closely related wall-lizard species as a model system. Molecular Phylogenetics and Evolution 125: 100–115. https://doi.org/10.1016/j.ympev.2018.03.021
  • QGIS.org (2023) QGIS Geographic Information System. Open Source Geospatial Foundation Project. Available from: http://qgis.org
  • Rambaut A, Drummond AJ, Xie D, Baele G, Suchard MA (2018) Posterior Summarization in Bayesian Phylogenetics Using Tracer 1.7. Systematic Biology 67(5): 901–904. https://doi.org/10.1093/sysbio/syy032
  • Recuero E, Buckley D, García-París M, Arntzen JW, Cogălniceanu D, Martínez-Solano I (2014) Evolutionary history of Ichthyosaura alpestris (Caudata, Salamandridae) inferred from the combined analysis of nuclear and mitochondrial markers. Molecular Phylogenetics and Evolution 81: 207–220. http://doi.org/10.1016/j.ympev.2014.09.014
  • Robin VV, Sinha A, Ramakrishnan U (2010) Ancient Geographical Gaps and Paleo-Climate Shape the Phylogeography of an Endemic Bird in the Sky Islands of Southern India. Gadagkar S (Ed.). PLoS ONE 5(10): e13321. https://doi.org/10.1371/journal.pone.0013321
  • Rohais S, Moretti I (2017) Structural and stratigraphic architecture of the Corinth Rift (Greece): An integrated onshore to offshore basin-scale synthesis. In: Lithosphere dynamics and sedimentary basins of the Arabian plate and surrounding areas. Springer, 89–120. Available from: https://doi.org/10.1007/978-3-319-44726-1_5
  • Ronquist F, Teslenko M, Van Der Mark P, Ayres DL, Darling A, Höhna S, Larget B, Liu L, Suchard MA, Huelsenbeck JP (2012) MrBayes 3.2: efficient Bayesian phylogenetic inference and model choice across a large model space. Systematic biology 61(3): 539–542. https://doi.org/10.1093/sysbio/sys029
  • Shackleton JD, Follows MJ, Thomas PJ, Omta AW (2023) The Mid-Pleistocene Transition: a delayed response to an increasing positive feedback? Climate Dynamics 60(11–12): 4083–4098. https://doi.org/10.1007/s00382-022-06544-2
  • Shen C-C, Miura I, Lin T-H, Toda M, Nguyen HN, Tseng H-Y, Lin S-M (2025) Exploring Mitonuclear Discordance: Ghost Introgression From an Ancient Extinction Lineage in the Odorrana swinhoana Complex. Molecular Ecology 34(10): e17763. https://doi.org/10.1111/mec.17763
  • Simon C, Buckley TR, Frati F, Stewart JB, Beckenbach AT (2006) Incorporating molecular evolution into phylogenetic analysis, and a new compilation of conserved polymerase chain reaction primers for animal mitochondrial DNA. Annual Review of Ecology, Evolution, and Systematics 37: 545–579. https://doi.org/10.1146/annurev.ecolsys.37.091305.110018
  • Simon C, Frati F, Beckenbach A, Crespi B, Liu H, Flook P (1994) Evolution, weighting, and phylogenetic utility of mitochondrial gene sequences and a compilation of conserved polymerase chain reaction primers. Annals of the entomological Society of America 87(6): 651–701. https://doi.org/10.1093/aesa/87.6.651
  • Song H, Mariño-Pérez R, Woller DA, Cigliano MM (2018) Evolution, Diversification, and Biogeography of Grasshoppers (Orthoptera: Acrididae). Insect Systematics and Diversity 2(4): 3. https://doi.org/10.1093/isd/ixy008
  • Song H, Amédégnato C, Cigliano MM, Desutter-Grandcolas L, Heads SW, Huang Y, Otte D, Whiting MF (2015) 300 million years of diversification: elucidating the patterns of orthopteran evolution based on comprehensive taxon and gene sampling. Cladistics 31(6): 621–651. https://doi.org/10.1111/cla.12116
  • Stefanidis A, Kougioumoutzis K, Zografou K, Fotiadis G, Tzortzakaki O, Willemse L, Kati V (2024) Mitigating the extinction risk of globally threatened and endemic mountainous Orthoptera species: Parnassiana parnassica and Oropodisma parnassica. Insect Conservation and Diversity 18(1): 54–68. https://doi.org/10.1111/icad.12784
  • Stefanidis A, Kougioumoutzis K, Zografou K, Fotiadis G, Willemse L, Tzortzakaki O, Kati V (2025) Distribution Patterns and Habitat Preferences of Five Globally Threatened and Endemic Montane Orthoptera (Parnassiana and Oropodisma). Ecologies 6(1): 5. https://doi.org/10.3390/ecologies6010005
  • Trifinopoulos J, Nguyen L-T, von Haeseler A, Minh BQ (2016) W-IQ-TREE: a fast online phylogenetic tool for maximum likelihood analysis. Nucleic Acids Research 44(W1): W232–W235. https://doi.org/10.1093/nar/gkw256
  • Ullrich B, Reinhold K, Niehuis O, Misof B (2010) Secondary structure and phylogenetic analysis of the internal transcribed spacers 1 and 2 of bush crickets (Orthoptera: Tettigoniidae: Barbitistini). Journal of Zoological Systematics and Evolutionary Research 48(3): 219–228. https://doi.org/10.1111/j.1439-0469.2009.00553.x
  • Uluar O, Yahyaoğlu Ö, Başıbüyük HH, Çıplak B (2023) Taxonomy of the rear-edge populations: the case of genus Anterastes (Orthoptera, Tettigoniidae). Organisms Diversity & Evolution 23(3): 555–575. https://doi.org/10.1007/s13127-023-00602-1
  • Weekers PH, De Jonckheere JF, Dumont HJ (2001) Phylogenetic relationships inferred from ribosomal ITS sequences and biogeographic patterns in representatives of the genus Calopteryx (Insecta: Odonata) of the West Mediterranean and adjacent West European zone. Molecular Phylogenetics and Evolution 20(1): 89–99. https://doi.org/10.1006/mpev.2001.0947
  • Wicker V, Ford M, Gawthorpe RL, Skourtsos E, Kranis HD, Kerouédan L, Muravchik M (2024) Transition from late-Miocene syn-orogenic extension to Plio-Pleistocene Corinth rifting in the southern Hellenides, northern Peloponnese, Greece. Tectonics 43(6): e2023TC007964. https://doi.org/10.1029/2023TC007964
  • Willemse F (1980) Three new species and some additional notes on Parnassiana Zeuner from Greece (Orthoptera, Ensifera, Decticinae). Entomologische Berichten 40(7): 103–112.
  • Willemse L, Kleukers R, Odé B (2018) The grasshoppers of Greece. EIS Kenniscentrum Insecten & Naturalis Biodiversity Center, Leiden, 439 pp.
  • Willemse L, Tilmans J, Kotitsa N, Trichas A, Heller K-G, Chobanov D, Odé B (2023) A review of Eupholidoptera (Orthoptera, Tettigoniidae) from Crete, Gavdos, Gavdopoula, and Andikithira. ZooKeys 1151: 67–158. https://doi.org/10.3897/zookeys.1151.97514
  • Xia X (2018) DAMBE7: New and improved tools for data analysis in molecular biology and evolution. Molecular Biology and Evolution 35(6): 1550–1552. https://doi.org/10.1093/molbev/msy073
  • Zhang J (2015) GMYC Web Server: web interface for single- and multi-threshold GMYC species delimitation. The Exelixis Lab. Available from: https://species.h-its.org/gmyc (September 16, 2025)

Supplementary material

Supplementary material 1 

Files S1–S13

Kotitsa N, Borissov SB, Chobanov DP (2026)

Data type: .zip

Explanation notes: File S1: Localities, specimens, collection data and Genbank accession numbers of sequences used in the phylogenetic analyses [.xlsx file]. — File S2: Substitution models and parameters used for Maximum Likelihood and Bayesian Inference analyses for the concatenated NAD2-COI-ITS matrix [.pdf file]. — File S3: Localities, specimens, xenocanto accession numbers, song recording and collection data used for the song analysis of Parnassiana [.xlsx file]. — File S4: Phylogenetic tree of the genus Parnassiana, inferred from the NAD2 genetic marker using Bayesian Inference analysis. Node labels indicate posterior probabilities [.pdf file]. — File S5: Phylogenetic tree of the genus Parnassiana, inferred from the COI genetic marker using Bayesian Inference analysis. Node labels indicate posterior probabilities [.pdf file]. — File S6: Phylogenetic tree of the genus Parnassiana, inferred from the ITS1–5.8S–ITS2 genetic marker using Bayesian Inference analysis. Node labels indicate posterior probabilities [.pdf file]. — File S7: Phylogenetic tree of the genus Parnassiana and a number of Platycleidini outgroups, inferred from the concatenated NAD2–COI–ITS matrix using Bayesian Inference analysis. Node labels indicate posterior probabilities [.pdf file]. — File S8: Phylogenetic tree of the genus Parnassiana, inferred from the concatenated NAD2–COI–ITS matrix using Maximum Likelihood analysis. Node labels indicate posterior probabilities [.pdf file]. — File S9: Time-calibrated tree of the genus Parnassiana and a number of Platycleidini outgroups, inferred from the concatenated NAD2–COI–ITS matrix. Node labels indicate mean node ages based on the NAD2 molecular rate for micropterous Orthopterans calculated by (Chang et al. 2020) (0.014256 subs/s/my). 95% HPD intervals and posterior probabilities for each node can be found in Supplementary file 10: BEAST table [.pdf file]. — File 10: Statistics (estimated age, 95% HPD interval min, 95% HPD interval max, posterior probability) of time-calibrated trees inferred from the concatenated NAD2–COI–ITS matrix. Three calibration schemes are included: 1. Calibration point at 4-3.5 mya at Node 3 to mark the separation of Peloponnese from the mainland; 2. NAD2 molecular rate for micropterous Orthopterans calculated by (Chang et al. 2020) (0.014256 subs/s/my); 3. NAD2 molecular rate for Pholidopterini calculated by (Çiplak et al. 2022) (0.018 subs/s/my). Node ID corresponds to the numbered nodes marked on the in-text timetree of Fig. 4 [.xlsx file]. — File 11: Includes qualitative and quantitative characteristics for the Parnassiana songs for each putative taxon, such as echeme length, microsyllable position and microsyllable number [.xlsx file]. — File 12: Includes the output of the GMYC species delimitation analysis for the NAD2 marker for Parnassiana as provided in the web server https://species.h-its.org [.pdf file]. — File 13: Includes the output visualisation and tables of the ABGD and ASAP species delimitation analyses for the NAD2 and COI markers for Parnassiana as provided in the SpartExplorer web platform https://spartexplorer.mnhn.fr [.pdf file].

This dataset is made available under the Open Database License (http://opendatacommons.org/licenses/odbl/1.0). The Open Database License (ODbL) is a license agreement intended to allow users to freely share, modify, and use this dataset while maintaining this same freedom for others, provided that the original source and author(s) are credited.
Download file (8.70 MB)
login to comment