To confirm the molecular identification of the five newly recorded populations (Escúllar, Felix, Torre de Maro, Río de la Miel, and Jete), we sequenced one to four specimens from each locality (Table 1). Mitochondrial gene fragments of the cytochrome oxidase subunit I (CoxI) and 16S ribosomal RNA (16S) were targeted for this purpose.

Genomic DNA was extracted from femoral muscular tissue using the Qiagen DNeasy extraction kit (Qiagen, Valencia, CA, USA). Polymerase chain reaction (PCR) was used to amplify CoxI and 16S using the set of primers LCO 1490/COI-H (Folmer et al. 1994; Machordom et al. 2003) and 16Sar/16Sbr (Palumbi et al. 1991), respectively. PCR conditions followed those described in Sánchez-Vialas et al. (2020). The PCR products were sequenced in both directions by Sanger Sequencing at Macrogen (Madrid, Spain).

Sequences were reviewed, assembled, and edited in Geneious v.11.0.18 and aligned using MAFFT (Katoh and Toh 2008). Phylogenetic analyses were conducted in BEAST v.1.10.4, following the dataset and methodology of Sánchez-Vialas et al. (2023), with the three newly generated sequences into the analysis incorporated.

Additionally, a phylogenetic network was generated using SplitsTree v.6.3.41 (Huson and Bryant 2006) based on combined CoxI and 16S mitochondrial sequences for all specimens possessing both markers. Although these markers reflect the maternal lineage, the network approach helps identify conflicting signals in the data – such as incomplete lineage sorting or historical gene flow – that may not be evident in traditional tree-based analyses. This method is particularly useful for detecting conflicting signals that might be obscured in traditional phylogenetic trees, offering a more comprehensive view of the evolutionary history.

We also considered the nuclear ITS2 allele data from Sánchez-Vialas et al. (2023) as additional evidence to support taxonomic decisions.

To investigate a possible effect of isolation by distance (IBD) within B. insignis sensu lato, we used the CoxI dataset to derive a genetic distance matrix (Table SS1) using MEGA X (Kumar et al. 2018). We selected CoxI because it is the most comprehensively represented gene in our matrix, shared by all sequenced specimens except one (N = 30). Genetic distance matrices were calculated using uncorrected genetic distances.

Geographic distances were calculated from georeferenced locality data (Table 1), obtained using Google Earth Pro. These distances were converted to Euclidean distances in kilometers using the distHaversine function from the geosphere package (Hijmans et al. 2017).

To assess the correlation between geographic and genetic distances, we conducted Mantel tests using the mantel function from the R package Vegan (Oksanen et al. 2013). Analyses were performed separately for each matrix, with results visualized through scatterplots.

A total of 182 specimens were examined, including 168 housed at the Museo Nacional de Ciencias Naturales (MNCN-CSIC, Madrid, Spain) and 14 at the Museo de Ciencias Naturales de Tenerife (MUNA, Sta. Cruz de Tenerife, Spain). Of these, 126 specimens were preserved in ethanol, while 56 were dry-preserved (Table S2). Additionally, detailed photographs of the holotype (by monotypy) of Meloe insignis Charpentier, 1818, deposited at the Museum für Naturkunde-Berlin (MNB, Berlin, Germany), were examined. These images were kindly provided by Bernd Jaeger, curator of Coleoptera at MNB. The holotype of Meloe insignis is labeled as follows: “Hist.-Coll. (Coleoptera), Nr.28477, Meloe insignis Charp., Zool. Mus. Berlin [with label, printed]; insignis, 1798. Typ. Charp* [white label, handwritten]; Type [red label, printed]; [QR code with former MFN collection URI].

Ethanol-preserved and dry-mounted specimens were examined under a stereomicroscope. Selected male specimens were rehydrated in water prior to genitalia extraction. For dry-mounted specimens, male genitalia were mounted on a cardboard using dimethylhydantoin formaldehyde (DMHF) resin and pinned adjacent the specimen. In ethanol-preserved specimens, male genitalia were stored in an Eppendorf tube filled with 96% ethanol, kept within the vial of the specimen. Measurements were taken from images captured using a Leica MZ16A stereomicroscope, fitted with a Leica DFC550 camera and processed with the software LAS v.4.3. Terminology for the male genitalia follows Selander (1966).

We conducted geometric morphometric analyses (GM) focusing on the shape variation of the last five antennomeres (VII–XI) in relation to four mitochondrial lineages of Berberomeloe insignis identified by Sánchez-Vialas et al. (2023): B1 (Murcia), B2 (eastern-central Almería), C1 (Granada), and C2 (south-central Almería). Male and female antennomeres were analyzed; antennomeres VII, IX, and XI showed sexual dimorphism, while antennomeres VIII and X did not. The analysis was based on 92 adult specimens, preserved dried or in ethanol. Among these, 55 were males, distributed across lineages B1 (13 specimens), B2 (15 specimens), C1 (13 specimens), and C2 (14 specimens), and 37 were females, distributed across lineages B1 (7 specimens), B2 (8 specimens), C1 (13 specimens), and C2 (9 specimens) (Table S3). Digital images of the inner side of left antennae (or right antenna, if left absent) were captured using a reflex camera coupled with a macro lens, mounted on a stable camera stand. For antennomeres VII–X, four anatomical landmarks were used, while for antennomere XI, seven anatomical landmarks were used for males and six for females (Fig. 1). Landmarks were digitized using the R package “StereoMorph” v.1.6.2. (Olsen and Westneat 2015) within R v.3.6.0 (R Core Team 2021). The coordinates (X, Y) obtained from landmark digitization were subjected to generalized Procrustes analyses (Gower 1975) to normalize for differences in position, orientation, and scale. Mean shapes were then computed for each lineage, and principal component analysis (PCA) was performed on these shapes. PCA was conducted on the coordinates projected into the linear tangent space using the function gm.prcomp. The first two principal components were used to visualize the variation in antennomere shape. Procrustes ANOVA with 1000 permutations was used to assess significant differences between clades.

Geographic records of B. insignis were compiled from previously published data (mainly García-París 1998; García-París et al. 1999, 2003; García-París and Ruiz 2008, 2011; Sánchez-Vialas et al. 2023), unpublished records from fieldwork (samplings carried out by the authors in the last five years), and observations taken from iNaturalist (www.inaturalist.org) (Table S4). From iNaturalist records, only those with high accuracy (N = 37) were used, excluding those with coordinate errors exceeding 1 km. Each iNaturalist observation was carefully assigned to its corresponding taxonomic unit based on cephalic coloration patterns.

The potential geographic distribution ranges of the five different lineages of Berberomeloe insignis were inferred using MaxEnt 3.4.1 (Phillips et al. 2017). MaxEnt estimates the interaction between environmental variables and species presence to create suitability models for a geographic area. Analyses were based on all known localities (Table S4). Nineteen bioclimatic layers were retrieved from the WorldClim 2.1 database (www.worldclim.org) at a 30 second resolution (Fick and Hijmans 2017). To avoid collinearity among predictors, we performed a Pearson correlation analysis in R using the package ‘corrgram’ (Wright 2006). When two variable showed a correlation coefficient higher than 0.8, only the one considered biologically more relevant for the species was retained. Based on this criterion, the following subset was selected: bio7 – Temperature Annual Range, bio8 – Mean Temperature of Wettest Quarter, bio11 – Mean Temperature of Coldest Quarter, bio14 – Precipitation of Driest Month, bio15 – Precipitation Seasonality (Coefficient of Variation), bio16 – Precipitation of Wettest Quarter. Analyses in MaxEnt were conducted using the cross-validation testing procedure recommended for small datasets (Philips et al. 2006). The MaxEnt output was set to Cloglog, which estimates the probability of presence between 0 and 1. The outputs were further edited with QGIS 3.22 Białowieża (www.qgis.org) to generate distribution maps.

Potential ecological interchangeability among the different lineages within the two main clades was evaluated by measuring niche overlap using Schoener’s D-metrics (Wooten and Gibbs 2012). Pairwise comparisons were conducted, and niche identity tests were performed using 100 randomized pseudoreplicates to determine if the Schoener’s D-values were statistically different (one-tailed test) than expected under the null distribution. Niche conservatism or divergence was inferred by comparing the niche overlap values with the null distribution of overlap values. Niche divergence is supported when the overlap values are smaller than the null distribution. These analyses were conducted with the R package “ENMTools” (Warren et al. 2021).

Distinguishing between interspecific and intraspecific variation is a challenging task (Pyron et al. 2020, 2024). To address species delimitation, we adopted the general lineage (unified) species concept (de Queiroz 1998, 2007), which considers a species as “separately evolving metapopulation lineages”. Based on this framework, we recognize species when populations show consistent patterns in both mitochondrial and nuclear DNA sequences (as demonstrated by the nuclear allelic network of ITS2 in Sánchez-Vialas et al. 2023) and exhibit concordant morphological traits. These multiple lines of evidence suggest that the lineages in question are evolving independently and have achieved reproductive isolation with minimal likelihood of future interbreeding, especially when the lineages involved are parapatric or sympatric.

The newly generated Cox1 (657 bp) and 16S (510 bp) sequences have been deposited in GenBank (Table 1). We recovered the same topology for Berberomeloe insignis sensu lato as in Sánchez-Vialas et al. (2023). Each of the newly sampled populations corresponds with its respective phenotypically diagnosable lineages (Fig. 2). Following the clade denominations of Sánchez-Vialas et al. (2023), the sequenced specimen from Escúllar (Sierra de los Filabres) is nested within the mitochondrial DNA (mtDNA) clade B2, the specimen from Felix (Sierra de Gádor) belongs to clade C2, and those from Río de la Miel, Torre de Maro, and Jete belong to clade C1.

The phylogenetic network generated using SplitsTree revealed a well-defined and cohesive lineage (C) encompassing populations from Granada, Málaga (C1), and western Almería (C2). Within this lineage, the specimens from western Almería showed a relatively genetic differentiation from western areas (C1). The clade C is clearly separated from the others by a long branch, indicating limited interaction with other lineages. In contrast, the remaining clades (B) exhibited a more intricate pattern, with high genetic differentiation but multiple reticulations (nexus points) within the network, suggesting extensive historical gene flow among them (Fig. 3).

The Mantel test results revealed low correlation between geographic and genetic distances across populations of Berberomeloe insignis sensu lato. Although the correlation was statistically significant (Mantel statistic r = 0.1165, P = 0.033), it remained notably weak (Fig. 4). This suggests that genetic differences among lineages are not primarily driven by geographic distance.

GM analyses of male antennomeres revealed a complex pattern of morphological variation (Figs 5, 6). Antennomere XI exhibited significant differences across most clades (Table S5), with the exception of the lineages of Central-Eastern Almería and Campo de Dalías (lineages B2 and C1 respectively, sensu Sánchez-Vialas et al. 2023), which showed a considerable overlap (Fig. 5). Antennomere X displayed substantial overlap among most clades, except for the one from Murcia (lineage B1 sensu Sánchez-Vialas et al. 2023), which showed significant differentiation. Antennomere IX followed the pattern of antennomere XI, with strong differentiation among lineages from Málaga-Granada, Murcia, and Central-eastern Almería, but with considerable overlap between the Central-Eastern Almería and Campo de Dalías clades. Antennomere VIII distinguished the Murcia lineage but did not differentiate among the other clades. Antennomere VII displayed significant differentiation between the Málaga-­Granada clade and other clades, with less pronounced differences between the Central-eastern Almería and Murcia clades (P = 0.02) and between the Central-eastern Almería and Campo de Dalías clades (P = 0.04) (Fig. 6).

In females, antennomeres showed substantial overlap among clades, except for the Murcia clade (B1), which exhibited clear differentiation from the others (Fig. 7). Notably, antennomere XI revealed a significant distinction between the C1 clade (Málaga-Granada) and the B2 clade (Central-Eastern Almería) (Table S4).

PCA for males indicated that PC1 and PC2 accounted for 45.37% and 14.87% of the variability for antennomere XI, 65.63% and 22.69% for antennomere X, 54.97% and 27.81% for antennomere IX, 67.65% and 19.75% for antennomere VIII, and 46.18% and 26.13% for antennomere VII. For females, PC1 and PC2 explained 31.9% and 21.9% of the variability for antennomere XI, 77.44% and 11.19% for antennomere X, 72.37% and 14.89% for antennomere IX, 66.16% and 18.66% for antennomere VIII, and 62.87% and 20.11% for antennomere VII.

In summary, all major clades (B1, B2, C1) can be distinguished based on variation in antennomeres, with the exception of the potentially hybrid populations represented by the mitochondrial subclade C2 (Campo de Dalías and southern Sierra de Gádor), which consistently exhibit significant overlap with other clades.

Head coloration patterns of the newly reported specimens are geographically congruent with their corresponding mitochondrial lineages (Sánchez-Vialas et al. 2023). Populations from Málaga, Granada, and the southwestern edge of Almería (clade C1 of Sánchez-Vialas et al. 2023) displayed two symmetrical red-orange blotches on the temples, slightly extending toward the frons (but not reaching the frons midline). Eastern Almería and Pasillo de Fiñana Valley (between Sierra Nevada and Sierra de los Filabres) populations (clade B2 of Sánchez-Vialas et al. 2023) exhibited three blotches: two symmetrical ones over the temples that do not extend towards the frons and a central V- or Y-shaped blotch on the frons. The Murcia phenotype (clade B1 of Sánchez-Vialas et al. 2023) resembles that from B2 but lack the central blotch on the frons. The remaining phenotype, located in Campo de Dalías and the southern slope of Sierra de Gádor (clade C2 of Sánchez-Vialas et al. 2023), resembles B1 but with wider separation of the temporal blotches at the occiput. Occasional, inconspicuous red blotches (very small and circular on both sides of the frontal midline, either isolated or connected to the temporal blotches) were also observed in some specimens from the southern areas of Sierra de Gádor.

The species distribution models (SDMs) performed well for the lineages of B. insignis, as indicated by high Area Under the Curve (AUC) values: Murcia (clade B1): AUC = 0.962; Eastern Almería and Pasillo de Fiñana Valley (clade B2): AUC = 0.950; Málaga, Granada, and the southwestern edge of Almería (clade C1): AUC = 0.983; Campo de Dalías and southern Sierra de Gádor (clade C2): AUC = 0.985. In addition, the SDM for all main lineages of B. insignis sensu lato yielded an AUC of 0.926.

The environmental variables that contributed most to the models varied by clade: for clade B1, bio16 (51.4%) and bio14 (15.9%); for clade B2, bio16 (50.2%) and bio15 (30.5%); for clade C1, bio6 (50.6%) and bio7 (25.3%); and for clade C2, bio14 (30.9%) and bio7 (30.4%). Across all lineages of B. insignis, bio16 (61.3%) and bio15 (14%) were the primary contributing variables.

Niche overlap, measured using Schoener’s D-metric, was low between clade C1 (Málaga, Granada, and the southwestern edge of Almería) and other clades, and relatively low to moderate among the different clades of B. insignis sensu lato (Table 2). All comparisons were statistically significant, indicating that the lineages are not ecologically interchangeable (Figs S1–S6).

Sánchez-Vialas et al. (2023) highlighted strong genetic isolation between the westernmost lineage (C1) and the other lineages, supported by nuclear and mitochondrial DNA. Considering this criterion, the morphological differences among the main units within B. insignis, and the species concept adopted (see material and methods), we propose recognizing two distinct species within the current conception of B. insignis sensu lato: one comprising lineage C1, which ranges from eastern limit of Málaga province (new faunistic records), southern Granada, and western limit of Almería province; and a second species comprising two subspecies distributed across Murcia (lineage B1) and eastern Almería (lineage B2).

According to this new taxonomic framework, since B. insignis sensu lato includes two candidate species, it is necessary to assign the original name to one of the two species. Charpentier (1818: 258) described Meloe insignis from Spain (type locality: “Spanien”) without precise locality. From the original description it can be concluded without any doubt that the species was described based on a single specimen, which constitutes the holotype fixed by monotypy (art. 73.1.3 of the International Code of Zoological Nomenclature; ICZN 1999). The type specimen illustrated in the original description shows a diagnostic head coloration pattern: two symmetric, reddish blotches on the temples, close to each other at the occiput, with no frontal blotch (Fig. 8A), confirmed by examining detailed photographs of the holotype (Fig. 8B, C). Thus, the phenotype of the name-bearing type (holotype) is characteristic and exclusive of the populations of Murcia province, corresponding to the B1 lineage. Therefore, we assign the name Berberomeloe insignis ( = Meloe insignis) to B1 lineage populations. As no synonyms exist for B. insignis, the new species (C1) and subspecies (B2) recognized here require, in turn, new names.

However, the situation becomes more complex for the populations of the lineage C2. Historical genetic introgression appears to explain the observed cyto-nuclear discordance in the populations from southwestern Almería (lineage C2). While mitochondrial DNA places these populations closer to the westernmost lineage (C1), nuclear ITS2 data aligns them more closely with the central-eastern Almería and Murcia populations (lineage B). Although geographically restricted to a small area in southern Almería, specimens from lineage C2 exhibit broad morphological variation, overlapping with their parapatric lineages. C2 populations are nevertheless distinguishable by a singular head coloration (see Sánchez-Vialas et al. 2023: figs 13, 14), yet some specimens from the southern slopes of Sierra de Gádor exhibit variability that may reflect a legacy of past hybridization. However, the available molecular data are limited, with only a few C2 individuals genotyped for mitochondrial DNA and ITS2, and consequently a significant portion of genetic diversity in these populations remains uncharacterized. Given these uncertainties, together with the absence of any fully diagnostic morphological trait, there is currently insufficient evidence to assess the extent of reproductive isolation between lineage C2 and its parapatric species (see discussion).

The current taxonomic catalogue for B. insignis should be modified with the following descriptions and diagnosis:

Over recent decades, the rise of molecular techniques has catalyzed a wave of cryptic species discoveries, leading to a deeper understanding of biodiversity and speciation processes (Witt et al. 2006; Adams et al. 2014; Pérez-Ponce de León and Poulin 2016; Pola et al. 2023). The recognition of cryptic diversity is increasingly significant across a range of taxonomic groups (Dufresnes et al. 2019; Sainz-Escudero et al. 2022; González-Miguéns et al. 2020). However, species initially identified as cryptic may reveal distinct morphological traits upon closer examination (Horsáková et al. 2019; Korshunova et al. 2017; Pola et al. 2023), suggesting reconsideration of their cryptic, or rather, pseudocryptic, status.

However, our study emphasizes the need to explore morphological diversification in species that may not be cryptic but are incompletely revised (Cuesta-Segura et al. 2023; Recuero and Caterino 2024). Our findings indicate that morphological and genetic changes observed across the studied populations are not part of clinal variation, as evidenced by distinct and abrupt morphological differences between them. Split tree analysis supports the previous hypothesis regarding their phylogenetic relationships, revealing that the easternmost lineage is well differentiated from the others. While moderate to strong isolation by distance (IBD) signals are typically expected in populations subjected to high gene flow, our results show only a slight overall trend toward IBD. The low Mantel statistic indicates that genetic differences among lineages are not primarily driven by geographic distance. This observation underscores the likelihood that other processes, such as incipient speciation, play a significant role in the diversification of this group, especially given the high divergence observed between clades. Based on the evidence, we identified and described B. nazarisp. nov. and B. insignis trisanguinatusssp. nov.

Our niche models show low overlap between B. nazari and other lineages, suggesting that it occupies a distinct climatic niche compared to the other taxa, including lineage C2 (Table 2; Figs S1, 18). The niche overlap among subspecies of B. insignis and the lineage C2 was relatively low to moderate, implying that while they share some climatic areas, notable ecological differences persist. Interestingly, despite phylogenetic and geographic proximity between B. insignis insignis and B. i. trisanguinatus, we did not observe specimens with intermediate characteristics that might indicate gene flow. This suggests that although these subspecies partially overlap in their climatic preferences, other factors, possibly ecological, behavioral, or genetic, may be maintaining their phenotypic differences near their contact zones.

In contrast, the absence of morphological differentiation and the presence of cyto-nuclear discordances in populations belonging to the lineage C2 highlight the need for further research to evaluate the extent of reproductive isolation among the involved lineages. At present, our data allow for outlining two alternative hypotheses that could explain the observed patterns of genetic and morphological discordances: (1) C2 populations would form a cohesive unit (sensu Coyne and Orr 2004), resulting from historical hybridization between lineages C and B, that led to an ancient mitochondrial capture, resulting in a hybrid species C2, with a closely related but differentiated C mtDNA lineage but with nuclear markers closely related to lineage B; Or, alternatively (2) C2 populations would not be cohesive, because they were the result of an old but still ongoing hybrid zone between lineage B2 and lineage C1, suggesting that B and C species coexisted in the region and that our sampling does not represent the extent of the ongoing hybridization. With our current sampling we cannot discriminate among these two alternatives. Further studies are needed to understand whether the processes described led already to speciation (hypothesis 1) or if we are facing a geographically unstable hybrid-zone (hypothesis 2). In the same geographic area, Andújar et al. (2012, 2014) presented a similar situation for the carabid beetle Carabus baguenai Breuning, 1926, a species resulting from ancient hybridization of nearly parapatric taxa, indicating that hybridization is not such a rare phenomenon in this geographic region. Since our data is incomplete we consider the C2 Berberomeloe problem still an open question for future research, and therefore, we refrain from making any taxonomic decisions regarding lineage C2. Additional sampling in zones between the ranges of B. insignis and B. nazari, combined with more robust genomic analyses (e.g., ABBA–BABA tests or F_ST estimates), will be necessary.

The addition of the newly described taxa in this study makes the province of Almería a focal point for Berberomeloe diversity (Fig. 18), including four taxa: B. indalo, B. insignis trisanguinatus, B. nazari, and B. tenebrosus (Sánchez-Vialas et al. 2020, 2023; this work). Most of these taxa exhibit parapatric distributions, with only B. indalo and B. insignis trisanguinatus found in syntopy at a few localities (García-París et al. 1999; this work). Even in a relatively small area such as the Sierra de Gádor, B. tenebrosus and the C2 lineage of B. insignis sensu lato, are parapatric, with the C2 lineage inhabiting lower altitudes (known populations found below 900 m a.s.l.) on southeastern slopes, while B. tenebrosus occupies higher altitudes (above 900 m, mostly between 1500 and 2000 m) on the western side (authors pers. obs.). A similar pattern is observed in Sierra de los Filabres, near Escúllar and El Haza de Riego, where B. tenebrosus and B. i. trisanguinatus occur less than 1 km apart. In this region, B. tenebrosus is confined to higher elevations and appears to be more common in its range (authors pers. obs.).

These findings present promising opportunities to study fine-scale ecological interactions across Berbero­meloe lineages, particularly the potential roles of competitive exclusion or altitude-based habitat partitioning, highlighting the need for further research on the mechanisms shaping species coexistence. In this regard, to refine their distribution boundaries, further sampling is needed between Campo de Dalías/Sierra de Gádor and La Rábita/Laujar de Andarax to the west (eastern limit of B. nazari and western limit of lineage C2, minimum distance: 29 km), and between Campo de Dalías/Sierra de Gádor and El Palmer/Aguadulce to the east (western limit of B. i. trisanguinatus and eastern limit of lineage C2, minimum distance: 12 km), to assess whether gradual phenotypic transitions occur. Similarly, additional sampling is needed for populations between the known localities of B. i. trisanguinatus and B. i. insignis (minimum distance: 18 km).

The observed variation in male antennal morphology across distinct lineages of Berberomeloe suggests a potential role for sexual selection and reproductive isolation in driving divergence (Carson and Bryant 1979; Kawano 2003). Secondary sexual traits often evolve due to sexual conflict or mate choice (Andersson 1994; Fricke et al. 2010), which can establish reproductive barriers and ultimately contribute to speciation processes (Gavrilets 2014).

We identified distinct morphological extremes in male antennomere XI. Males of B. nazari have a wider antennomere, notched in its apical area (except some specimens of the C2 lineage), while B. i. insignis and B. i. trisanguinatus exhibit slender antennomeres, typically unnotched. Other antennomeres also show notable differences: X and VIII are significantly elongated in B. i. insignis compared to other taxa, while IX differs significantly across B. i. insignis, B. i. trisanguinatus, and B. nazari. Additionally, B. nazari has a shorter antennomere VII compared to the other taxa. Specimens from the C2 lineage exhibit a remarkable morphological variability in antennomeres XI and IX, along with head coloration, with antennal morphology overlapping that of B. i. trisanguinatus and B. nazari. A plausible explanation for the increased morphological variability in specimens from the C2 lineage is the existence of past or ongoing gene flow (Fuzessy et al. 2014; Enciso-Romero et al. 2017; Grant and Grant 2019). Notably, males and females of C2 show comparable variation, whereas in established species male variation is typically more constrained, likely due to sexual conflict or mate choice. This pattern suggests that C2 could represent a hybrid zone rather than a cohesive species and highlighting an avenue for future genomic investigation.

Male antennomeres V, VII, IX, and XI show a strongly dentate apical border on the inner side, covered by very short setae, distinct from the setae on the rest of the antennomere (authors pers. obs.). This specialized setation might play roles in searching for females, mate recognition, and courtship behaviors. In Berberomeloe, courtship begins with antennation, where males use their antennae to scrutinize the last female’s abdominal tergites through vibratory movements (Bologna 1988). The variation in male antennae across these taxa may reflect adaptations related to mating behaviors, mate recognition, or environmental interactions (Greenfield 2002; Jayaweera and Barry 2017). However, non-adaptive processes, such as genetic drift, cannot be excluded. Field and experimental studies are needed to explore the behavioral and ecological significance of this antennal morphological diversity.

Female antennal morphology remains largely conserved, except in B. i. insignis, which shows marked differentiation in antennomeres VII, IX, and X. This elongation of antennomeres in females of B. i. insignis is particularly intriguing, especially given the apparent lack of selection pressures for antennal function in females, in contrast to the strong selection males face during mating.

The challenge of identifying and naming evolutionary units, whether intraspecific or at the species level, is crucial for raising awareness and ensuring the protection of distinct populations (Liu et al. 2022). Each unique lineage may harbor evolutionary traits that contribute to broader biodiversity, underscoring the need for precise taxonomic recognition. In this context, some populations might be distinct enough genetically or ecologically to warrant protection but do not yet meet the usual standards for species designation under commonly applied species concepts (i.e., due to the absence of clear reproductive barriers). Addressing this taxonomic challenge requires consideration of subspecific classification frameworks (Dufresnes et al. 2023, 2024). Berberomeloe presents a complex taxonomic situation due to its significant genetic and morphological diversity. The geographic proximity of the main units of Berberomeloe, along with the existence of narrow unsampled regions between them, suggests a predominantly parapatric distribution rather than clear allopatry. While parapatric taxa are often expected to show gradual phenotypic transitions due to gene flow (Mayr 1963), Berberomeloe instead exhibits sharp phenotypic distinction between units, suggesting potential barriers to gene flow or strong selective pressures maintaining these differences.

Conservation strategies should not only focus on taxa but also on protecting regions where key evolutionary processes – such as hybridization and local adaptation – drive and maintain biodiversity. Hybrid zones and other areas of ongoing diversification are crucial for preserving the mechanisms that generate evolutionary novelty and ecological resilience. This is particularly urgent in rapidly deteriorating landscapes such as southeastern Spain (Sánchez-Vialas et al. 2023). In this regard, Sierra de Gádor and its surroundings emerge as a critical focal area for conserving evolutionary processes in Berberomeloe.

Our findings emphasize the urgent need for intensified taxonomic exploration, even for taxa previously considered well-documented. Regions like Campo de Dalías and Sierra de Gádor, which are under severe threat from greenhouse expansion, highlight the critical importance of precise taxonomic identification to inform effective conservation efforts (Liu et al. 2022; Dufresnes et al. 2023; Recuero et al. 2023). Beyond their significance as centers of endemism, these areas also serve as key sites for understanding evolutionary processes, such as historical and ongoing hybridization. Preserving both the species and the evolutionary dynamics shaping them is crucial, reinforcing the need to integrate these factors into conservation planning.