Abstract

Comprising nearly 500 species worldwide, Cheilosia Meigen is the largest genus of Syrphidae in the Palaearctic region. Within Cheilosia, phenotypic diversity has been assessed in different species groups, including the group of Cheilosia longula (Zetterstedt, 1838). However, the phenotypic variability of Cheilosia soror (Zetterstedt, 1843), a highly variable member of the C. longula group, has never been assessed in western Europe. In the present work, morphological and molecular analyses were conducted to assess the phenotypic variability found in 300+ specimens of C. soror from the Iberian Peninsula. A total of 16 variable characters were identified and defined for the C. soror morphology, with the highest variation found in the colour of the mesonotum pilosity and the metatibia colour. Morphological variation was assessed against molecular variation based on two molecular markers, one mitochondrial, the 5 end of the cytochrome c oxidase subunit I (COI-5), and one nuclear, the large subunit ribosomal DNA (28S). Phylogenetic analyses rendered trees with topologies in disagreement with the defined morphological variation. Two haplotypes were identified amongst the analysed specimens of C. soror, together with a haplotypic variant exclusive to the Iberian region. Potential distributions were used to identify unexplored areas of occurrence of C. soror and other species of the C. longula group in the Iberian Peninsula. Unassessed areas of occurrence of C. soror should be surveyed in the future to confirm the absence of hidden taxonomic diversity within the range of phenotypic variation for this species. Phenotypic variation of the other two Iberian species of the C. longula group, C. longula and C. scutellata (Fallén, 1817), was also assessed to find that they are species with less-variable morphology than C. soror and with molecular characters in accordance with other conspecific populations in Europe. New distributional data are provided for C. soror and C. scutellata from Spain, and a leg abnormality is identified for the first time in C. soror.

1. Introduction

The genus Cheilosia Meigen, 1822 (Syrphidae: Eristalinae) comprises black to brown syrphids with shiny thorax and abdomen covered with black, yellow, or red hairs [1]. It is one of the largest syrphid genera, with over 475 species described worldwide [2, 3] and the largest in the Palaearctic Region [4, 5]. Although saprophagy is present among Cheilosia larvae [6, 7], most species have phytophagous larvae, tunnelling in stems and roots of plants, leaf-mining, and feeding on cambium and sap [6]. Cheilosia larvae can even be associated with deciduous trees [7].

The genus Cheilosia is taxonomically complex due to its high species richness [8], the high level of cryptic speciation [1], and the low degree of interspecific variation reported in many cases [9]. Thus, to facilitate Cheilosia classification and identification, authors have often divided the genus into different morphological species groups (e.g., [10, 11]), sometimes supported by molecular characters too (e.g., [12, 13]). For example, the group of Cheilosia longula (Zetterstedt, 1838), which is the target of the present study, was defined by Barkalov [14] based on the male genitalia and redefined by Clauβen and Ståhls [12]. The last mentioned provided the most updated key to the European species of this group to date. The species belonging to the C. longula group share bare eyes, connected antennal sockets (nonseparated by the lunula), bare face (excluding microtrichia), and tibiae pale at least basally and apically. This group comprises the following species worldwide: Cheilosia flavissima Becker, 1894, C. longula, Cheilosia pallipes (Loew, 1863), Cheilosia scutellata (Fallén, 1817), Cheilosia soror (Zetterstedt, 1843), and Cheilosia thessala Clauβen & Ståhls, 2007. Larvae of the C. longula group are known to feed in the fruiting bodies of fungi (Basidiomycetes) [15].

An integrative approach using different sorts of data (e.g., morphological, molecular, and morphometric) has proven useful in the resolution of complex taxonomic problems in the Syrphidae (e.g., [1622]). Many of these taxonomic problems originate from observing phenotypic variation in a certain species. For example, by using morphology, morphometry, and DNA characters of Chrysotoxum Meigen, 1803, Nedeljković et al. [23] found that, behind the variation of the colour of the hairs of the thoracic scutum in Chrysotoxum vernale Loew, 1841, there were two separate species, C. vernale with yellow hairs and Chrysotoxum montanum Nedeljković & Vujić with black hairs. However, the high phenotypic diversity observed in Cheilosia vernalis (Fallén, 1817) by Ståhls et al. [1] was found to be not in accordance with the variability of the 3 end of the mitochondrial gene cytochrome oxidase subunit I (COI-3), and only genetic haplotypes were defined within this species.

Also within Cheilosia, Finnish populations of Cheilosia naruska Haarto & Kerppola, 2007 showed both phenotypic diversity and low genetic variability in a study by Milankov et al. [24]. Wing geometric morphometrics was used to assess patterns of phenotypic differentiation in the group of Cheilosia variabilis (Panzer, 1798) and Cheilosia melanopa redi Vujić, 1996 was upgraded to the species level, as Cheilosia redi [25]. Within the C. longula group, there is only one published work about an analysis of phenotypic diversity [26]. In the mentioned work, DNA, proteins, and geometric morphometry were used to assess the phenotypic diversity observed between two populations of an undescribed species, Cheilosia aff. longula. Although interpopulation divergence in wing shape and five haplotypes in total were identified for the analysed Finnish populations of C. aff. longula, the detected genetic variation was low and taxonomic diversity could not be identified within C. aff. longula [26]. This low genetic variation has great importance in the conservation of threatened and endangered species, considering that it could be correlated with lower evolutionary potential, as well as short-term decreases in fitness and long-term lack of adaptive flexibility, for example, to cope with climatic change [24, 26]. Regarding C. soror, two genetic haplotypes have been identified along with data on geometric morphometrics of wing in the only assessment known of the genetic/phenotypic diversity of (Balkan) C. soror [27].

As it is not straightforward to monitor syrphids due to their size or because it is often not possible to identify them in situ, potential distribution studies might well serve for assessing the conservation status of a species and its population dynamics or for characterising its most suitable habitat. These studies also play a key role in establishing new protected areas as a fundamental strategy for preserving biodiversity against environmental or anthropic pressures [28, 29]. In addition, potential distribution studies can be useful to unravel unexplored areas of occurrence of a species where hidden taxonomic diversity within the same species could be present.

With 57 species, Cheilosia is the genus with the most species recorded in the Iberian-Balearic region [30, 31]. Three species belong to the C. longula group: C. longula, C. scutellata, and C. soror. Preliminary examination of material of C. soror collected in different Spanish localities of the Iberian Peninsula revealed high morphological variability in different characters. However, it is known that phenotypic diversity within some Cheilosia species is not always congruent with taxonomic diversity (e.g., [1]). Our observations prompted us to formally assess this phenotypic variability and explore a hypothetic congruence between phenotypic and taxonomic diversities. Thus, the specific objectives of the present study were (1) to provide an account of the variable morphological characters observed in the studied Iberian material of C. soror; (2) to assess the genetic diversity of this material in correlation with the observed morphological variability; (3) to understand how other species of the C. longula group found in the Iberian Peninsula vary morphologically and genetically in comparison with C. soror; and (4) to define expected areas of occurrence of C. soror in the Iberian Peninsula with unexplored populations of this species suitable for future phenotypic assessments.

2. Material and Methods

2.1. Target Species: Distribution and Breeding Sites

The three species of the C. longula group found in the Iberian Peninsula were the targets of the present study, with a special focus on the highly variable C. soror. Species were identified with Clauβen and Ståhls [12] and Vujić [32]. Cheilosia longula occurs in most European countries, from Fennoscandia to the Pyrenees and from Ireland to Siberia. In Europe, the flight period of C. longula is from late June to October, with a peak in September [8]. Its larva has been found in fungal fruiting bodies of different species of Boletus L., Leccinum Gray, and Suillus Gray [6]. Cheilosia scutellata is widespread in the Palaearctic region, occurring from Fennoscandia to the Iberian Peninsula and North Africa and from Ireland to the Pacific coast. The flight period of C. scutellata goes from May to October, up to November, in southern regions of Europe [8, 33]. Its larva can be frequently found in the fruiting bodies of Boletus or Suillus, among many others [6, 33]. Cheilosia soror occurs from Fennoscandia south to North Africa and southern England to Siberia and Japan. The adults of C. soror fly from May to September, with a peak in June/July [8]. Its larva remains undescribed but is reported as having been found in truffles (e.g., Tuber), Leccinum, Suillus, or Xerocomus Quélet [15, 33, 34].

2.2. Examined Material

In total, 518 specimens of the C. longula group were analysed: C. longula: 10♂♂ and 2♀♀; C. scutellata: 98♂♂ and 89♀♀; and C. soror: 173♂♂ and 147♀♀. The examined material (see Supplementary material 1) includes specimens collected with hand net in different areas of mainland Spain, from Sierra Nevada (south) to Picos de Europa (north), between 2019 and 2021, but also unpublished Iberian specimens deposited in the Spanish collections “Colección Entomológica de la Universidad de Alicante” (CEUA-CIBIO) and “Museo Nacional de Ciencias Naturales, Madrid” (MNCN). A few specimens from eastern Europe were also examined for comparison from the entomological collection of the Faculty of Sciences, University of Novi Sad, Serbia (FSUNS). To avoid unnecessary repetitions in the examined material, the collection acronym is mentioned only when the specimen/s is/are deposited in a collection other than the CEUA-CIBIO. Some specimens were collected with Malaise traps and are indicated with “T.Ma.” in the examined material list (Supplementary material 1). All specimens from the CEUA-CIBIO bear a barcode label with a unique identifier code.

2.3. Morphological Study

The morphology of the involved species was studied in detail, including the male genitalia. Morphological terminology follows Thompson [35] except for “hairs” that are used instead of “pile.” For the male genitalia, terminology follows Clauβen [36]. Specimens were observed using a binocular stereomicroscope Leica© M80. Photographs were taken by Focus Stacking with a Leica© DFC 450 camera attached to a Leica© M205 C binocular microscope, using Leica© Application Suite (LAS X) software. Male genitalia were prepared following the protocol in Ricarte et al. [37] and stored in microvials with glycerol. Male genitalia were drawn on a Stone© TH light tablet from stacks of photos produced with the equipment mentioned above and software. As a result of the morphological assessment, a table with the list of characters found to be variable in C. soror was produced, with an indication of the sex and a definition of the range of variability for each character. Regarding C. longula, due to the low number of specimens available from the Iberian Peninsula, it has not been possible to conduct a proper study of morphological variability.

2.4. Molecular Markers

A total of 25 specimens (Supplementary material 2) were used for sequencing the cytochrome c oxidase subunit I (COI-5), and 13 of those for sequencing the nuclear large subunit ribosomal DNA (28S). DNA was extracted from both the right mesoleg and metaleg of single-pinned individuals using the NZY gDNA isolation kit following manufacturer’s protocol. PCR amplification of the 5 region of the COI-5 gene was performed with the universal primers LCO1490 (5-GGTCAACAAATCATAAAGATATTGG-3) and HCO2198 (5-TAAACTTCAGGGTGACCAAAAAATCA-3) [38]. PCR amplification of the 28S was performed with the universal primers F2 (5-AGAGAGAGTTCAAGAGTACGTG-3) and 3DR (5-TAGTTCACCATCTTTCGGGTC-3) [39]. PCR products were visualized with an electrophoresis process in a 1% agarose gel and sequenced at Macrogen Inc. (Macrogen, Spain).

For the COI-5 and 28S, PCR reactions were performed in 25 μl reactions containing 1 μl DNA extract, 0.5 μl of both respective primers (0.2 μM), 0.2 μl of DNA polymerase (1 unit) (0.4 μl for 28S sequencing), 1 μl MgCl2 (2 mM), 2.5 μl 1× buffer, 1 μl of dNTPS (0.4 mM), and ultrapure water. COI-5 and 28S amplifications followed Ståhls and Barkalov [13] and Belshaw and Quicke [40], respectively. The annealing temperature was modified in 28S amplification, varying between 48 and 50°C.

2.5. Phylogenetic Analyses

In total, 25 COI-5 and 13 28S sequences from Iberian specimens were obtained successfully. All species of the C. longula group, except Cheilosia thessala (no sequences available online) and Cheilosia pallipes (North American species), were compared by an analysis of sequences. All information of sequences obtained in this study is available in GenBank (see Supplementary material 2). Phylogenetic analyses were carried out under maximum likelihood (ML) and Bayesian inference (BI) frameworks. Sequence alignment was performed manually with the program AliView v. 1.27 [41] for ML analysis and automatically with MAFFT online service [42] and checked with AliView for BI analysis. The final matrix had a length of 657 bp (ML) and 1215 bp (BI). A maximum likelihood tree for the COI-5 gene was built in MEGA7 v.7.0.26 [43] with 1000 bootstrap replications, using a general time reversible (GTR) model with gamma distribution (+G) which was the best evolutionary model proposed by the corrected Akaike information criterion (AICc, 3275.323) using the software jModelTest v.2.1.10 [44]. Ferdinandea cuprea (Scopoli, 1763) was used as an outgroup from GenBank (http://www.ncbi.nlm.nih.gov/genbank).

Pairwise distance analyses between the studied specimens were conducted in MEGA7, and the computations were run by default using the p-distance model with 1000 bootstrap replicates (see Supplementary material 3). A Bayesian consensus tree from a concatenated matrix of COI-5 and 28S sequences from Iberian specimens of C. longula, C. scutellata, and C. soror was performed in MrBayes v.3.2.7 [45]. We ran the Bayesian inference analysis with four chains for 1 million of generations, sampling every 100 generations, in order to diminish standard deviation under 0.01. With a view to select the best-fit nucleotide substitution model for each partition of the concatenated matrix, we employed PartitionFinder2 v.2.1.1 [46] under the Bayesian information criterion (BIC) and using the greedy algorithm [47]. Two subsets dividing the sequences of both markers were established in order to set the best-fit nucleotide substitution models, but both partitions obtained the same best-fit model, HKY; therefore, we decided to work with only one partition containing both genes. The 50% majority rule consensus tree was rooted with Ferdinandea aurea Rondani, 1844 sequences obtained from GenBank and then edited in FigTree v.1.4.4 (http://tree.bio.ed.ac.uk/software/figtree).

2.6. Potential Distribution Analysis

Maps of the actual distribution of C. longula, C. scutellata, and C. soror in the Iberian Peninsula (plus Balearic Islands and some North African localities) were produced with the software QGis [48] based on all the specimens included in Supplementary material 1, as well as data on bibliographic (nonrevised) records following Ricarte and Marcos-García [30] and Van Eck [49, 50]. Species records without a geographical coordinate given in the specimen labels and/or in literature were assigned the most precise coordinate possible by using Google Earth. All these data were also used for building the predictive distribution models, made by using the free-access software MaxEnt [51]. Predictive distribution models for C. scutellata and C. soror and his possible larval hosts (fungi) were also provided and displayed in a combined format in which resulting likelihood equates to the summation of both models’ likelihood. The hosts were selected from literature amongst the putative and confirmed ones having fruiting bodies during the maximum activity peak of the mapped Cheilosia species. The chosen host fungi for one or other of the two studied Cheilosia species were Boletus, Gyroporus Quélet, Lactarius Persoon, Leccinum, Pholiota Fries, Suillus, Tuber P. Micheli ex F. H. Wiggers, and Xerocomus. Geographical records of the mentioned fungi were downloaded from http://www.gbif.org.

For building the initial models, 19 bioclimatic variables were used describing current climate obtained from WorldClim dataset [52, 53] in 30 arc sec resolution (Table 1). To reduce collinearity, the modelling was performed in two stages. In the first part, all bioclimatic variables were used for model building, while in the second step, modelling was conducted using only the strongest predictor variables, i.e., those showing more than a 10% of contribution in the initial model [28, 29]. Those variables listed below had a percentage contribution higher than 10% in the first step and therefore were used for the debugged model (d. m.). Same procedure was applied for every mentioned fungus genus. For the combined larva-host models, each included host genus is indicated with their +10%-contribution variables from its own d. m.: C. scutellata d. m. (; ): C. scutellata larva-host d. m. ((Boletus: ); (Gyroporus: ; ); (Lactarius: ); (Leccinum: ); (Pholiota: ); (Suillus: ; ); (Tuber: ; ; ; and ); (Xerocomus: ; ; and )). C. soror d. m. (; ): C. soror larva-host d. m. ((Boletus: ); (Tuber: ; ; ; and )).

Resulting potential distribution models were combined and displayed together with the actual distribution maps for C. scutellata and C. soror. For C. longula, only the map of the actual distribution was presented (see Section 3.3). For reference of location of Spanish provinces, see Ricarte and Marcos-García [30].

3. Results

3.1. Phenotypic Variation Assessment

Of the three studied species, C. soror showed the highest morphological variability. The list of variable characters and the description of their range of variation are reported in Table 2. Male genitalia also displayed morphological variability in certain parts. The typical male genitalia of C. soror is illustrated in Figures 1(a)1(d). One of the 77 dissected genitalia of C. soror showed an atypical morphology (CEUA_S16) (Figures 1(e)1(h)). Still, any other external difference was found in comparison with the other studied males of C. soror. The outlier genitalia had the ventral lobe of the gonostylus shorter and thicker apically, the surstylus was both more slender and round-shaped (with more angular edges in typical genitalia), and the average surstylus lamella comprised almost the entire surstylus length (up to 70% of the surstylus length in typical genitalia). Within the morphological assessment of C. soror, a female collected in Sierra de Aitana (Alicante, Spain) in 2021 was found to have anomalous metabasitarsi, which were strongly swollen (CEUA_S69) (Figures 2(a) and 2(b)). The maximum width of the metabasitarsi was twice the apical width of the metatibia.

The few studied individuals of C. longula showed well-defined characters without any variation, not even in those featured as the most variable traits of C. soror (Table 2). For C. scutellata, the most remarkable variable characters were the colour of the mesonotum hairs in both sexes and the colour of all tibiae, mainly in females. Both traits showed similar variation as for C. soror, although in the case of the mesonotum hairs, the variation was not as extreme. The studied genitalia of C. longula and C. scutellata showed no intraspecific variability (Figure 3).

3.2. Phenotypic Variability vs. Genetic Variability in Cheilosia soror

Cheilosia soror was the most variable species in morphology (see Phenotypic Variation Assessment). So, this morphological diversity was tested against genetic diversity based on the markers COI-5 (mitochondrial gene) and 28S (nuclear gene). A set of specimens with various variable characters and from different geographical origins within Spain were analysed (Table 3).

Compared with a mitochondrial DNA sequence of Drosophila yakuba Burla, 1954 [54], we obtained a 657 bp fragment of COI spanning nucleotide positions 1516 to 2173. Two main haplotypes were registered (A, B), along with 4 variants for haplotype A (A1a, A1b, A2a, and A2b) (Figure 4), for all the used C. soror sequences (Supplementary material 2).

Within haplotype A, the A1a variant was found to be exclusive to the Iberian Peninsula. Specimens with the A1a variant (CEUA_S16; S17; S18; S19; S21; S24; S25; S69; S112; S115; S118; S119; and S121) showed a dark pattern both in the mesonotum hairs and tarsi, but a few specimens were lighter in these characters (Table 3). A1b variant was represented by a single specimen from Germany (BIOUG17467-E05), but we did not examine this specimen’s morphology, same for the specimen with the A2a variant (BIOUG16275-D03). The A2b variant belonged to a single female from Spain and was exclusive to “Picos de Europa” (León) (CEUA_S114) and also exclusive from the Iberian Peninsula. However, this female from “Picos de Europa” had no outstanding morphological traits (Table 3). Regarding haplotype B (CEUA_S20; S23; S110; S111; J_Skevington_22751; BIOUG16618-D06; and BIOUG17375-E03), it occurs both in Spain and central and eastern Europe with specimens, at least the Iberian ones, showing a pale pattern both in the mesonotum hairs and tarsi, opposite to the A1a variant (Table 3).

3.3. Potential Species Distributions

The actual occurrence points for each of the three Iberian species of the C. longula group are represented in Figure 1. Based on a low number of Iberian records, C. longula was excluded from the potential distribution analyses. The projected distributions of C. scutellata (Figure 5(a)) and C. soror (Figure 5(b)) were roughly similar. Cheilosia scutellata had an average probability of occurrence in the Pyrenees above 80%; meanwhile, C. soror showed a probability near 60-70%. On the other hand, C. soror had an average probability of occurrence in the southeast of the Iberian Peninsula (provinces of Alicante and Valencia) near 50%, while the average probability of C. scutellata remained around 30-40%. According to the maps, the northern half of the Iberian Peninsula has greater options to host populations of C. scutellata with an average likelihood of 70%, as well as northeastern Mediterranean shores and southern parts of Spain with an average likelihood of 40-50%. Regarding C. soror, the northern half of the Iberian Peninsula has also shown high probabilities of occurrence, averaging 60% of likelihood, but with greater affinity for Mediterranean environments than C. scutellata (50% vs. 30% average probability of occurrence, respectively).

4. Discussion

4.1. Molecular vs. Morphological Characters

Cheilosia soror is a morphologically highly variable species (Tables 2 and 3). One of the most variable characteristics of adult morphology is the relation between the presence of yellow and black colours on femora and tibiae, as well as the presence of black and yellow hairs on mesonotum, abdomen, and legs. After examination of all the listed specimens (Table 2; Supplementary material 1) we can conclude that the only constant characters in C. soror are, in fact, the diagnostic characters of this species, i.e., the presence of hairs on the anterior anepisternum in both sexes and the bright orange enlarged basoflagellomere of females, as well as the bare eyes, face with indistinct hairs, at least the base and the apex of tibiae pale, and lunula not dividing antennal sockets [32]. Variation of morphological characters was not following the species phenology as specimens captured on similar dates showed completely different morphologies (Table 3). Our results suggested that there is no correlation between the phenotypic variation and environmental factors, as it has been demonstrated for the syrphid genus Eupeodes Osten Sacken, 1877 [55].

We concluded that the atypical genitalia of a C. soror male caught in Sierra de Aitana (CEUA_S16) (Figures 1(e)1(h)) is a result of extreme intraspecific variability. The specimen shared the haplotype A1a with other specimens with typical genitalia (Figures 1(a)1(d)). No other morphological character of the specimen CEUA_S16 showed phenotypic singularity, and the idea of this specimen representing a separate taxon was discarded. Variable genitalia have been reported from different species in Syrphidae genera such as Heringia Rondani, 1856 or Merodon Meigen, 1803 [17, 56], and also regarded as intraspecific variability. The usually slight changes reported here for the surstylus of a typical C. soror genitalia (Table 2and Figures 1(a)1(d)) appeared to be heavier in the atypical one (Figures 1(e)1(h)), suggesting that it represents an extreme intraspecific variability or even a teratological feature.

Concerning the malformed female of Sierra de Aitana (CEUA_S69) (Figures 2(a) and 2(b)), the only known report of the same aberration in the C. longula group was provided by Ståhls and Claussen [57] as occurring in females of C. thessala and C. longula. This is the first time that a leg aberration is notified in C. soror. In addition, Clements [58] described a rare leg structure occurring in a female of Cheilosia albitarsis (Meigen, 1822). Although this reported deformation in C. albitarsis does not resemble the abnormality of our female, Clements cited some examples of similar aberrations in the genera, including a personal comment of Steven Falk referring to “a highland population of the species C. longula throwing up females with swollen rear metatarsi.” The malformed C. soror female (CEUA_S69) was caught at an altitude of 1130 m in Sierra de Aitana, the highest mountainous formation in the province of Alicante.

Both maximum likelihood (Figure 6) and Bayesian consensus trees (Figure 7) showed no differences between analysed individuals of each species, well separated in different tree branches, but displaying a complete polytomy inside every branch. We hypothesized that the specimens found somewhat further apart from the polytomy were not different species since this divergence might be due to the change in a single nucleotide between their haplotypes. However, this single change was neither supported by the analysis of the pairwise distance (Supplementary material 3), which also corroborated the raised hypothesis with a divergence range of 0%-0.51% between C. soror specimens (Supplementary material 3), 0%-0.3% for C. longula, and 0%-0.6% for C. scutellata specimens.

Based on the evidence given by the molecular analyses, despite using both mitochondrial and nuclear markers, we could not find consistent genetic patterns in C. soror to support taxonomic differences in correlation with the observed morphological variation. Although there is a slight colour tendency in the two main haplotypes (A, “dark”; B, “pale”), we cannot confirm this variation as a general pattern due to the small number of individuals with defined haplotypes and the existence of colour exceptions in specimens of both haplotypes. Haplotype A variants also need to sample more individuals to ensure their taxonomic situation, especially the A1a and A2b variants, exclusive to the Iberian Peninsula.

The variety of haplotypes and respective variants in the Iberian Peninsula as well as the fact that some of them are shared with northern and eastern Europe might be caused by a speciation process due to Pleistocene glaciations and the existence of glacial refugia during this period in the Iberian region [59]. Pleistocene speciation has also been reported in some hoverfly genera, such as Eumerus and Merodon [60, 61], focused in the Balkan region, which is also recognised as Pleistocene glacial refugia. Furthermore, the Iberian exclusivity of haplotypes A1a and A2b may indicate that there is ongoing speciation in the Iberian Peninsula, but in an early stage, not allowing a clear molecular/morphological differentiation between taxa yet.

It might also be that the morphological changes observed in the studied specimens of C. soror are only intraspecific variability. Our inferred trees employed genetic sequences of individuals of different parts of Europe together with our sequences of specimens with different morphological traits from several Iberian locations. The results confirmed that the species concepts of C. scutellata and C. soror are the same all over Europe, with insufficient genetic differentiation. The low genetic variability revealed in the Iberian C. soror may be related with a lower evolutionary potential, as well as a decrease in fitness in the short term or lack of adaptive flexibility in the long run [26], common characteristics in threatened and endangered species, unlike the species studied in the present paper that are regarded as least concern in Europe/EU27 [6264].

Contrary to expectations due to its low genetic variation (Figures 57), the great morphological variability of C. soror and its widespread distribution across most of the Palaearctic region might be a result of phenotypic plasticity. Brought on by genetic and/or environmental factors, phenotypic plasticity increases fitness. It may aid dispersal and colonization [65, 66], explaining the wide distribution of C. soror and its local abundance and, somehow, reversing the negative consequences of low genetic variability. Further studies of how environmental factors influence the phenotypic plasticity of C. soror are desirable, as well as molecular studies implicating gene markers or haplotypes, to confirm the hypotheses put forward in this work.

4.2. The Implication of Potential Distributions in Forthcoming Studies

Although the sampling effort of the present study has been rather extensive, our predictive models showed that there are still many areas and locations in Spain with a high probability of harbouring species of the C. longula group that have not been sampled yet (Figures 1(a) and 1(b)), especially the northern provinces of Spain with high occurrence probabilities as La Rioja, Galicia, or Navarra. The areas previously indicated in the results as not prospected might be of great importance to continue assessing phenotypic and genetic diversities of C. soror, as well as the rest of the species of the C. longula group. An individual of C. scutellata was recently captured for the first time in Sierra Nevada National Park (Granada) [67] (Ballester-Torres et al. unpublished data), coinciding with our predictive model, which considered the locality as suitable for the species (40-70% likelihood) (Figure 1(a)), and therefore, demonstrating the model effectiveness. After three exhaustive surveys in the national park in different times of the year, only one individual was caught. This could explain the 40-70% likelihood rates of the predictive model in this area instead of percentages close to 100. The degree of coincidence between predictive models and reality could increase if we use other variables that consider the connection of hoverflies with flowering plants and their preference for open or closed areas, like plant type and its coverage percentage or land use [68, 69]. As an example, fungus genera employed in our potential distribution analysis have been proved to be very effective to make the model stronger and more accurate, as the expected distribution of the studied species of Cheilosia largely coincided with those of the fungi studied. These hoverflies appear closely related to fungi, as they are rarely found in locations where mushrooms are missing. Therefore, the survival of the species of C. longula group may depend largely on the development of fruiting bodies of the studied fungi, as suggested in the literature [6, 8, 15]. Years with low rainfall and humidity might drive to a decrease in the number of developing fruiting bodies and, therefore, in a decline in the number of new adults of the three species of the C. longula group [70]. Thus, the presence of the studied fungal hosts may be an almost certain indicator of the presence of the species of the C. longula group and vice versa. However, we have not always captured individuals in those places marked as optimal by the potential distribution models. For example, Lleida and La Rioja are considered as suitable for C. scutellata and C. soror by the predictive model (70-100% likelihood) (Figures 1(a) and 1(b)), but these species have never been recorded from these provinces [30]. The reasons why these species have not been recorded from both provinces can be diverse, for example, due to the absence of syrphid surveys at the appropriate time of the year, absence of fungal hosts (highly improbable due to the humid climate of these geographic areas), or just insufficient sampling effort in these areas.

All this being said, new sampling efforts may be worthwhile in coniferous and deciduous forests (e.g., Sierra de Francia, Salamanca; Selva de Oza, Huesca); deep valleys (e. g., Valle del Jerte, Cáceres; Ordesa y Monte Perdido, Huesca), high-altitude ecosystems (e.g., nonstudied locations of the Pyrenees as Valle de Ansó, Huesca, or the Navarrese Pyrenees; Sierra de Gredos), or riverside locations included in these potential areas, as they might be harbouring isolated populations resulting in new haplotypes or even different species.

Data Availability

The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

Ethical Approval

All applicable international, national, and international guidelines for the care and use of animals were followed. This article does not contain any studies with human participants performed by any of the authors.

Disclosure

The present work is part of Iván Ballester-Torres PhD thesis.

Conflicts of Interest

The authors declare that they have no conflict of interests.

Acknowledgments

This research was funded by the “Fauna Ibérica” project (PGC2018-095851-A-C65) of the Ministry of Science, Innovation, and Universities, Spain. Antonio Ricarte’s position (Ref. UATALENTO17-08) at the University of Alicante is funded by the “Vicerrectorado de Investigación y Transferencia del Conocimiento.” The positions of Zorica Nedeljković and Iván Ballester-Torres at the University of Alicante are funded by the above-mentioned Fauna Ibérica project. Thanks are due to Mercedes París for facilitating the access to the MNCN collection and arranging the loan of specimens. We are very grateful to Pablo Aguado-Aranda for helping with the molecular part of this study. We also thank Patrik Rada for his helpful comments and advice on GIS software.

Supplementary Materials

Supplementary 1. Supplementary material 1: list of examined material for each species of the C. longula group in the Iberian Peninsula.

Supplementary 2. Supplementary material 2: GenBank and/or BOLD accession numbers for the gene sequences used in this study. “CEUA” sequences are new.

Supplementary 3. Supplementary material 3: pairwise distance analyses between the studied specimens from maximum likelihood tree inferred from partial COI-5 gene sequences.