Chromosome preparations were obtained from 10 blind mole rat individuals from five localities of the S. microphthalmus sensu lato range (Fig. 1; Tables 1, S1). The sample comprised five individuals from four localities of the S. microphthalmus (2n = 60) range and five individuals (2n = 62, Spalax sp.) from the vicinity of Kislovodsk town (Stavropol Krai, Russia) in the North Caucasus. The standard technique for preparing chromosome slides from bone marrow was applied (Ford and Hamerton 1956). Routine staining of chromosomes was performed with 2% Giemsa. Differential G-banding was conducted using the method of Seabright (1971), and C-banding was carried out according to Sumner (1972). The stained chromosome slides were observed and photographed using a KEYENCE BZ-9000 (Keyence Corporation, Japan) microscope with immersion. Karyograms were assembled in Adobe Photoshop v.19.1.2. We followed the chromosome arrangement scheme proposed by L’apunova et al. (1974).

Synaptonemal complexes (SCs) were prepared by the microspreading technique following Matveevsky et al. (2021). For immunocytological analyses, chromosome spreads were mounted on poly-L-lysine–coated slides. The primary antibody used was rabbit anti-SYCP3 (1:500; Abcam, ab15093), which labels the axial and lateral elements of the SCs. The secondary antibody was Alexa Fluor 488–conjugated goat anti-rabbit IgG (Jackson ImmunoResearch). Immunostaining was performed according to standard procedures described by Matveevsky et al. (2021). Fluorescent signals were examined using an Axio Imager D1 microscope (Carl Zeiss, Jena, Germany). Measurements of SC lengths were performed using MicroMeasure software (Colorado State University, CO, USA).

Tissue samples (kidney, liver, or skin) for the study of the genetic variability of blind mole rats were obtained from the Collection of wildlife tissues for genetic research of the Koltzov Institute of Developmental Biology of the Russian Academy of Sciences (CWT IDB), state registration number 3579666. A total of 20 individuals of S. microphthalmus sensu lato from 11 localities were studied, including six from the North Caucasus (Fig. 1; Tables 1, S1). This sample included two previously karyotyped individuals from the Kursk region with 2n = 60 (Puzachenko and Baklushinskaya 1997; see our Table SS1). In addition, two specimens of the lesser blind mole rat Nannospalax leucodon (Nordmann, 1840) from Odessa Oblast, Ukraine (CWT IDB collection numbers 25248, 25249) were examined for comparative analysis.

Genomic DNA was extracted by the standard salt method (Aljanabi and Martinez 1997). The full-length cytochrome b (cyt b) mitochondrial gene and the interphotoreceptor retinoid binding protein (IRBP) nuclear gene were used as phylogenetic markers. Twenty individuals of S. microphthalmus sensu lato were sequenced for cyt b; of these, five individuals with 2n = 60 and six with 2n = 62 were also sequenced for IRBP (Table SS1). A polymerase chain reaction (PCR) was performed using specific primers SpalaxCBfw2 and SpalaxCBrv2 for cyt b (Németh et al. 2020), and F11 and R22_cric (Lebedev et al. 2018) for IRBP in 15 μL volume contaning 0.3 U HS-Fuzz polymerase and 6 μL 2.5X reaction buffer (Dialat Ltd., Moscow), 0.3 mM dNTP (Evrogen, Russia), 1 pM of each primer, 7.85 μL ddH₂O, and 1 μL DNA 30 ng/μL in a Veriti Thermal Cycler platelet amplifier (Applied Biosystems, Waltham, MA, USA). PCR was performed under the following conditions: initial denaturation for 3 min at 95 °C; 35 cycles of denaturation for 30 sec at 95 °C, annealing for 40 sec at 60 °C for cyt b or 50 sec at 60 °C for IRBP, extension 1 min 30 sec or 50 sec at 72 °C, and final extension for 7 min at 72 °C. Automated sequencing was performed using a NovaDye Terminator Cycle Sequencing Kit 3.1 (GeneQuest, Moscow) with the AB 3500 Genetic Analyzer (Applied Biosystems, Waltham, MA, USA), and the Nanofor-05 Genetic Analyzer (Syntol, Moscow) at the Core Centrum of IDB RAS. Sequence chromatograms were reviewed and manually edited using the SeqMan section of the Lasergene 11 package (DNASTAR, Madison, WI).

All newly generated DNA haplotype sequences have been deposited in the NCBI GenBank database under accession numbers PV012717–PV012730, PX613617, and PX613618.

For the phylogenetic analysis, the obtained datasets of cyt b and IRBP sequences were supplemented with three cyt b sequences and two IRBP sequences of S. microphthalmus available from GenBank (Table SS2). Sequences of the following species, retrieved from GenBank, were also used in phylogenetic analyses: six other Spalax species (Table 1) – S. antiquus Méhely, 1909 (3 cyt b and 2 IRBP); S. arenarius Reshetnik, 1939 (3 and 2); S. giganteus (2 and 1); S. graecus Nehring, 1898 (1 and 1); S. uralensis Tiflov & Usov, 1939 (2 and 1); and S. zemni Erxleben, 1777 (= S. polonicus Méhely, 1909) (5 and 2), as well as three Nannospalax species as an outgroup – N. ehrenbergi (Nehring, 1898) (3 and 1); N. leucodon (2 and 1); and N. xanthodon (Nordmann, 1840) (3 cyt b; Table SS2). Sequence alignments were performed with MUSCLE algorithm (Edgar 2004) in MEGA 12 software package (Kumar et al. 2024) and subsequently corrected manually.

Maximum likelihood (ML) phylogenetic analysis was performed in IQ-TREE v3.0.1 (Wong et al. 2025), while Bayesian inference analysis (BI) was performed in MrBayes v3.2.7 (Ronquist et al. 2012). Two separate analyses were conducted for cyt b: one based on a full-length (1140 bp) cyt b alignment, which excluded the shorter (582 bp) S. uralensis sequences, and another based on a reduced-length (582 bp) alignment that included sequences from all species. The best model for nucleotide sequence evolution under the Bayesian information criterion (BIC) was selected separately for the 1st and 2nd codon positions (partition 1) and the 3rd codon position (partition 2) of both cyt b alignments, as well as for IRBP alignment using ModelFinder (Kalyaanamoorthy et al. 2017). The following models were used to construct ML / BI trees: TPM2u+F+I / GTR+I for partition 1 of both cyt b alignments; TIM2+F+I+G4 / GTR+I+G4 for partition 2 full-length cyt b; TN+F+G4 / GTR+I+G4 for partition 2 of short cyt b alignment; and HKY+F / HKY for both partitions of IRBP. BI analysis was based on 2 million generations with each 5000^th retained in two MCMC chains with default other settings. The first 25% of the generated trees were discarded as burn-in. Node support for ML trees was assessed using ultrafast bootstrap (UFBoot; 10,000 replicates) and SH‑like approximate likelihood ratio test (SH‑aLRT). Nodes with SH‑aLRT ≥80% and UFBoot ≥90% were considered supported; SH‑aLRT of 70–80% and UFBoot of 80–90% were interpreted as marginal. Posterior probability (PP) values were used to estimate branch support for BI trees, with PP ≥0.90 considered significant. Trees were visualized using FigTree v1.4.3. (http://tree.bio.ed.ac.uk). Genetic differences were assessed by pairwise distances (p-distance) and Kimura 2-parameter (K2p) for full-length cyt b sequences in MEGA 12. Interspecific genetic distances were calculated for samples larger than 3 individuals.

A Median Joining (MJ) evolutionary network of full-length cyt b and IRBP haplotypes were constructed in the HaplowebMaker software (https://eeg-ebe.github.io/HaplowebMaker [accessed on 27 January 2025]; Spӧri and Flot 2020). All SpalaxIRBP sequences generated in this study or retrieved from GenBank were homozygous (containing no ambiguity codes), with the exception of S. giganteus, which possessed an ‘R’ at position 119. The latter sequence was phased manually.

The genetic variability full-length cyt b indices, including the number of haplotypes (H), haplotype diversity (Hd), and nucleotide diversity (π), were estimated for the S. microphthalmus (2n = 60) sample in DnaSP v.6.12.03 (Rozas et al. 2017). We also calculated Fu and Li’s F*, Tajima’s D and Fu’s Fs (Tajima 1989; Fu 1997) in DnaSP to assess historical population growth, decline, or stability. The expansion coefficient (S/k) was calculated to assess the differences between recent and historical population sizes, as the ratio of the number of variable sites (S) to the average number of pairwise nucleotide differences (k) (Peck and Congdon 2004).

A total of 142 intact skulls of adult blind mole rats belonging to six species of the genus Spalax were analyzed from 57 localities across the species ranges (Fig. 1; Tables 1, S1). Among these, the number of skulls belonging to S. microphthalmus from 25 localities and Spalax sp. (2n = 62) from the 5 localities were 43 and 18, respectively (Fig. 1; Table 1). A population sample of 122 skulls of adult S. microphthalmus (51 females, 60 males) from the Streletskaya Steppe (Central Tsernozemny Biosphere Nature Reserve, Kursk district, Kursk oblast, Russia) was used as additional material. Adulthood was determined by scoring the morphological features of skull structure, such as the development of crests and the chewing surface structure of molars (Topachevsky 1969). Topachevsky (1969: 71) suggested that sexual size dimorphism in mole rats is not pronounced and can only be detected in series-based material through the mean values of certain cranial traits. We preliminarily estimated the proportion of variance in cranial measurements attributable to sexual size dimorphism across all studied Spalax species (Table SS3a). On average, this contribution amounted to 6.55 ± 0.08%, which was more than an order of magnitude lower than the corresponding contribution of interspecific variation (68.8 ± 3.1%). Therefore, we combined males and females into a single sample.

The primary skull collections included in this study are housed in the Zoological Institute of the Russian Academy of Sciences (ZIN, St. Petersburg) and the Zoological Museum of Lomonosov Moscow State University (ZMMU, Moscow). The sample also included skulls of specimens previously collected in the region of interest (Kabardino-Balkaria), which are stored in the collections of Kabardino-Balkarian State University (KBSU, Nalchik) and the Tembotov Institute of Ecology of Mountain Territories of the Russian Academy of Sciences (TI, Nalchik). The total collection comprised skulls from both previously karyotyped animals (KBSU) and non-karyotyped specimens (TI, KBSU; Table SS1).

A preliminary comparison of the latter with the karyotyped specimens revealed similarities in diagnostic morphological features (e.g., the structure of the hard palate and the small size of adult individuals). At this stage of the study, this grouping was treated as a hypothesis, to be further tested using morphometric analysis.

Twenty-eight measurements were taken of each skull using digital callipers with an accuracy of 0.1 mm. The scheme of measurements used is shown in Figure 2. We used the ratio of a particular measurement to the maximum skull length in a number of cases to characterize relative sizes. In addition, the length of incisive foramen and the length and width of upper and lower molars were estimated in 9 specimens of S. microphthalmus and 6 specimens of Spalax sp. (2n = 62). The width and length of the chewing surface of a sample of 23 Spalax teeth from the Middle Pleistocene (MIS 9, 11: 337–300, 424–374 ka BP) paleosols of the Otkaznoye locality (44°18'56.97"N, 43°50'44.29"E, Stavropol Krai, Russia (Markova 2006; Bolikhovskaya et al. 2016)) were also measured: 3 M1, 4 M2, 4 M3, 6 m1, 4 m2, and 2 m3. We used the dental terminology of upper and lower molars according to (Topachevsky 1969; Sarika and Sen 2003; López-Antoñanzas 2012; Erten 2018), which are shown on molars of S. microphthalmus juveniles in Figure 3.

Common univariate statistical methods used included descriptive statistics (median, mean, standard error of the mean, and 95% confidence limits for the mean and percentiles), Pearson’s correlation, the Welch test for comparisons of group means under unequal variances (Welch 1947), and the nonparametric Mann–Whitney U test (Mann and Whitney 1947) as an alternative to the independent two-sample t-test when normality could not be assumed or when sample sizes were small (N ≈ 20 or less). We also applied two normality tests (Shapiro–Wilk and Lilliefors with Monte Carlo estimation) to our data, which did not reveal significant deviations from the null hypothesis in the vast majority of cases (SM2). As a precaution, in all cases we routinely ignored “statistical significance” at p > 0.01.

Two morphospaces (Puzachenko 2023) were developed using the nonmetric multidimensional scaling (NMDS; Kruskal 1964; Davison 1983) based on the matrices of Euclidean distances between all pairs of skulls (general skull size model, SZM) and Kendall’s (1938) tau-b (corrected for ties) associations (general skull shape model, SHM). The structure of the SZM and SHM models was studied using correlation analysis. For the species samples, the medians of the coordinates of the SZM and SHM models were calculated.

In the text, the NMDS axes were denoted as E1, E2, etc. for SZM models, and as K1, K2, etc. for SHM models. We used components of variance analysis (Solomon 2005) to estimate the a priori taxonomy effect on the model axis variations. The species sample centroid coordinates were used as variables in a cluster analysis.

Note, NMDS has several advantages over parametric dimensionality-reduction methods such as principal component analysis (PCA) or factor analysis (Minchin 1987; James and McCulloch 1990; Legendre and Legendre 1998). It does not rely on assumptions of normality, linear relationships among variables, or homogeneity of variances, which are often violated in morphometric datasets. NMDS operates on ranks of dissimilarities and can therefore accommodate any distance measure, making it more flexible for heterogeneous morphological characters. By preserving the relative ordering of inter-object distances rather than their absolute values, NMDS is better suited to capturing nonlinear patterns in the data and is less sensitive to differences in scale among variables. In addition, NMDS often performs more robustly with limited sample sizes and provides an ordination space that directly reflects morphological dissimilarity among specimens. However, we applied principal component analysis (PCA) to construct a “metamodel,” in which the variables were NMDS axes rather than skull measurements, in order to account for the correlations between the axes of the SZM and SHM models (i.e., correlations between variation in skull size and shape).

Standard statistics (mean, standard errors, etc.), normality tests, the Welch test, the Mann–Whitney U test, and PCA were carried out using PAST v. 3.12 (Hammer and Harper 2007). NMDS analyses were performed using NCSS v. 12 Statistical Software (Professional License [Perpetual]; ncss.com/software/ncss). Dendroscope ver. 3.8.3 (Huson and Scornavacca 2012) was used for the preparation of some dendrograms.

Among the ten examined S. microphthalmus sensu lato individuals, five possessed chromosome sets with 2n = 60, NFa = 116, and NF = 120 (Fig. 4). These specimens originated from the terra typica of the species in the Novokhopersk steppes (Topachevsky 1969), as well as from the Kursk region, the southern bank of the Don River delta, and the northern bank of the Manych-Gudilo Channel (Table SS1). This karyotype was similar to that previously described for S. microphthalmus (L’apunova et al. 1974). The other five blind mole rats were from the North Caucasus and had 2n = 62, NFa = 120, and NF = 124 (Fig. 5A). This chromosome set consisted of five pairs of medium-sized and small metacentrics, 12 pairs of progressively smaller submetacentrics, and 13 pairs of subtelocentrics, including two pairs of the largest elements of the set. The X chromosome was a large metacentric or submetacentric, while the Y chromosome was a small subtelocentric. This karyotype corresponds to that previously described from the Northern Caucasus (Dzuev and Shogenov 2004; Dzuev et al. 2025a), except for the ratio of metacentrics to submetacentrics. This discrepancy may be explained by ambiguity in the morphology of small bi-armed chromosomes with different degrees of spiralization.

Heterochromatin blocks were detected in the centromeric regions of all chromosomes in the 62-chromosome karyotype, except for the smallest submetacentric pair (Fig. 5C). The X chromosome had a large subcentromeric block, and the Y chromosome was almost entirely heterochromatic. This result was consistent with the textual description of C-banding provided by Dzuev and Shogenov (2004): “… pericentromeric heterochromatin blocks have been detected in nearly all chromosomes …”. The near-centromeric heterochromatin blocks were detected on all chromosomes in the normal karyotype of S. microphthalmus (2n = 60) and were consistent with the previously published C-banded karyotype (Puzachenko and Baklushinskaya 1997).

As expected, spermatocytes of the 60-chromosome mole rat (specimen #27269) exhibited 29 autosomal bivalents, allowing identification of 29 fully formed SCs and the sex bivalent at the pachytene stage (Fig. 6A). In the 62-chromosome mole rat (specimen #27268), 30 autosomal bivalents were formed, and accordingly 30 complete SCs and the sex bivalent were distinguished at pachytene (Fig. 6B).

In both individuals, the sex (XY) chromosomes were positioned at the periphery of the meiotic nucleus, forming a characteristic sex body (Fig. 6). Notably, during the transition from late zygotene to early pachytene, the SYCP3 axes of the sex chromosomes underwent reorganization, transforming into rounded, sometimes reticulated structures of variable size. These paired structures corresponded to the X and Y chromosomes. Such “lumpy” configurations were consistently observed in both male mole rats (Fig. 6).

In the total sample of 23 complete cyt b sequences (1140 bp) from S. microphthalmus sensu lato, thirteen distinct haplotypes were identified in S. microphthalmus (2n = 60) and two haplotypes in Spalax sp. (2n = 62). The 60- and 62-chromosome forms shared no common haplotypes. The mean nucleotide composition was as follows: A = 31%, T = 32%, G = 13%, and C = 24%. The dataset contained 106 variable positions, 96 of which were parsimony informative. A total of 107 mutations were identified, of which 18 were non-synonymous substitutions. Within the cytochrome b protein sequence, which contains 379 amino acids, 11 were found to be variable in the total sample of S. microphthalmus sensu lato. Six of these variable amino acids distinguished the 60- and 62-chromosome forms. Analysis of the nine obtained IRBP sequences (859 bp) revealed two haplotypes: one shared by four S. microphthalmus individuals (2n = 60), and the other by five Spalax sp. individuals (2n = 62; Table SS1). All studied individuals were homozygous for the IRBP gene.

The topologies of the ML and BI phylogenetic trees of both mitochondrial and nuclear genes were found to be similar (Figs 7, S1). While all studied Spalax species formed distinct, statistically well-supported branches on all trees, statistical supports for deep nodes in the ML trees were low or marginal (Fig. S1). The phylogenetic position of several clusters also remained unresolved in some BI trees. Specifically, the placement of S. giganteus and the combined S. giganteus / S. uralensis cluster were unsupported in all analyses. Furthermore, the positions of S. graecus and S. antiquus on the nDNA tree also lacked statistical support (Figs 7B, S1B). The well-differentiated branches of Spalax sp. and S. microphthalmus grouped together into a common clade in all analyses, albeit with low statistical support on short alignment mtDNA trees (Figs 7A, S1A).

The MJ network of cytochrome b (cyt b) haplotypes clusters the haplotypes of different species separately (Fig. 8). Species-specific haplogroups are separated by a relatively high number of mutational steps, ranging from 56 between S. arenarius and S. zemni to 88 between S. microphthalmus and Spalax sp. (2n = 62). Haplotypes of Spalax sp. are linked to the rest of the network via haplotype c13 (specimen 27257), collected at the southern margin of the S. microphthalmus range (Table SS1). Sequence analysis revealed that the connection between c13 and the Spalax sp. haplotypes was determined by two sites: an Adenine at position 495 and a Guanine at position 813. Both represent synonymous transitions within the S. microphthalmus cyt b sequence located at the third codon position and are specific to Spalax sp. This contrasted with other haplotypes of this S. microphthalmus lineage (c10–c12) found in southern populations.

The IRBP haplotype network also resolved haplotypes of different species into separate branches (Fig. 8B). The number of mutational steps separating S. microphthalmus and Spalax sp. was comparable to that separating S. giganteus and S. uralensis, as well as S. arenarius and S. zemni.

The average cyt b gene-based genetic distances between the blind mole rat species studied are presented in Table 2. All Spalax species exhibited significant genetic differentiation. The genetic distances between Spalax sp. and S. microphthalmus reached notional interspecific values (>5%) for mammals (Baker and Bradley 2006), comparable to the differences observed among other Spalax species, including both of these species and S. giganteus. Despite its broad distribution, S. microphthalmus exhibited remarkably low intraspecific genetic variability. Its cyt b gene-based genetic distances were comparable to those of closely related species with far more restricted ranges, such as S. arenarius, S. graecus, and S. antiquus (Rusin et al. 2024). This result was consistent with the previously reported low level of S. microphthalmus genetic variability of the control region of mtDNA (Matveeva et al. 2019). The genetic variability indices for the S. microphthalmus cyt b sample (N = 17) were: H – 13; H_d – 0.963; π – 0.0028; S – 20; k – 3.191. The values of neutrality tests were: Fu and Li’s F* = –2.39673, P>0.05; Tajima’s D = –1.82274, P<0.05 and Fu’s Fs = –7.878, P<0.01. The expansion coefficient (S/k) was 6.268.

The SZM and SHM models comprised four and three coordinates respectively (Table SS3a). The relative components of variance in coordinates E1 and K1 due to the conventional taxonomy in the genus Spalax were 83.7% and 80.4% respectively. The next most important coordinates that show a “taxonomic signal” were K2, E4, and K3. The proportion of explained variance attributed to the E1–E4 coordinates varied from 78% (DDL) to 97% (MSL); on average, it was 90%. Thus, the model satisfactorily describes variability in skull size. The high mean proportion of variance (0.81) of the skull measurements attributed to coordinate E1 indicates the dominance of skull variability associated with the overall skull size.

Key measurements of size variation were MSL, RSMIH, RSMEW, PL, MSW, FGMAW, INCW, MNDL, and others (Figs 9A, S2A, S3A; Table SS3b). Interorbital width (YW) correlated with both E2 (r = 0.67) and E3 (r = 0.65) coordinates. The projections of the species sample centroids onto coordinates E1 and E2, arranged along E1 according to the disparities in species-specific overall skull size. In this sequence, the skull of Spalax sp. (2n = 62) occupies the outermost position as the smallest on average. At the opposite end of the row were the largest members of the genus Spalax: Spalax uralensis and S. giganteus. Because the sample of S. arenarius consisted mainly of female skulls, in Figure 9A, the centroid of this species is located in close proximity to the centroid of the small Spalax sp. (2n = 62).

Coordinate K1 of the SHM showed a weak but significant correlation with MSL (r = 0.56). This indicates the presence of an allometric pattern, i.e., the shape of a skull depends partly on its general size. K1 most strongly correlates with indices TUL/MSL (r = 0.80), TDL/MSL (r = 0.73), INCW/MSL (r = 0.79), and RSMEW/MSL (r = 0.77), i.e., with the relative sizes of tooth rows, the width of the upper incisors, and rostrum width (Figs 9C, S2B, S3C; Table SS3). Coordinate K2 demonstrated the strongest correlation with indices RSMIW/MSL (r = 0.61), RSMIH/MSL (r = 0.54), MNDH/MSL (r = 0.60), CNDW/MSL (r = 0.58), and YW/MSL (r = –0.57).

Cluster analysis of sample centroids over all coordinates of the SZM or SHM showed a strong separation of S. uralensis and S. giganteus (Puzachenko 1993) from all the other species of the genus by size and shape of the skull (Fig. 9B and Fig. 9D, respectively). This finding is consistent with the placement of these species on a distinct branch of the phylogenetic tree (Fig. 7A, B). Spalax sp., as the smallest representative of the genus, occupied a near basal position in the classification tree for all other species (Fig. 9B). It was closest to S. graecus, then to S. microphthalmus, and was well differentiated from both S. zemni and S. arenarius on the tree reproducing similarity in overall skull shape between the taxa (Fig. 9D).

At last, we applied PCA of the coordinates of the SZM and SHM together to combine information from both models (Figs 9E, S4; Table SS4). In accordance with this synthesis, Spalax sp. shares a “subclade” together with S. graecus and S. microphthalmus. We note, however, that bootstrap support for differentiation was generally low, except for the division into two groups (S. uralensis + S. giganteus vs. all other species), which can be explained by substantial morphological overlap among the species (Fig. S3).

In general, the results of the study of craniometric variability agree with the findings of genetic studies and do not contradict the hypothesis of radial speciation of modern Spalax from a common ancestor. Furthermore, the observed convergence of Spalax sp. with S. microphthalmus and S. graecus is consistent with the evolutionary cyt b haplotype networks of Spalax (Fig. 8).

The North Caucasian blind mole rat is a typical member of the genus Spalax, based on a combination of the following craniometric characteristics (Topachevsky 1969): the supracondylar foramen is always absent; the pharyngeal tubercles (tuberculum pharingeum laterale) are narrow, elongated as narrow ridges on the edges of the basisphenoid bone (corpus ossis sphenoidalis) and the basilar part of the occipital bone; the basisphenoid bone and the basilar part of the occipital bone are approximately level with each other and the skull base is not “fractured”; the corono-alveolar notch (incisura corono-alveolaris) is strongly developed and sharp; the fossa between the corono-alveolar and corono-condylar notches is deep; and the angular process of the mandible is attached to the alveolar process (see SM2 for more detail).

Spalax sp. differs from all other extant species of Spalax primarily in its smaller skull size. Nevertheless, there is no morphometric hiatus between most skull measurements of this form and those of species (which is typical for Spalax species) other than very large S. uralensis and S. giganteus (Table S5). Statistically significant differences (p = 0.0021 to < 0.0001) between Spalax sp. and all other species were detected for MSL, ZW, CPSH, INCW, TUL, TDL, and FGL.

Variation in several indices (VAR/MSL, %) demonstrates specific characteristics of Spalax sp. (Fig. S5, Table S5), including a relatively high rostrum, the widest skull base (mastoid width) and interorbital breadth, the longest mandible, and a relatively short glenoid fossa. Thus, although the absolute maximum rostrum height in Spalax sp. is comparable to that of other species (12.9 ± 0.21 mm, 11.26–14.37 mm), its relative height (25.6 ± 0.35%, 22.7–27.9%) is significantly greater (p = 0.0056 to < 0.0001) than in all other species except S. uralensis. The interorbital width in Spalax sp. (7.7 ± 0.12 mm, 6.65–8.81 mm) does not differ substantially from this measurement even in very large species such as S. giganteus or S. uralensis, in which it ranges from 6.2 to 10.0 mm. Only S. graecus is characterized by a wider interorbital region on average (8.2 ± 0.14 mm, 7.4–9.2 mm). However, the relative index of this measurement in Spalax sp. (15.4 ± 0.28%, 13.1–17.2%) is significantly higher (p = 0.0052 to < 0.0001) than in all other species, except S. graecus (14.9 ± 0.31%, 13.3–17.6%). An even clearer example of non-allometric variation is provided by mastoid width, which reflects the overall dimensions of the skull base. In Spalax sp., the absolute mastoid width (26.3 ± 0.23 mm, 24.37–28.05 mm) is significantly smaller (p = 0.004 to < 0.0001) than in all species except S. arenarius. However, the relative mastoid width reaches 52.2 ± 0.32% (50.3–55.1%) of the maximum skull length, which is significantly greater (p = 0.0044 to < 0.0001) than in the remaining Spalax species.

According to our data, Spalax sp. has a considerably shorter foramen incisivum compared with S. microphthalmus (Table 3), and possibly compared with other Spalax species, which remains to be verified in future studies. The average length of the IFL is 18.7% of the length of the upper diastema in Spalax sp. and 24.2% in S. microphthalmus. It is surprising that the relative length of the IFL in Spalax sp. corresponds to values typical for species of the genus Nannospalax – N. ehrenbergi (Coşkun et al. 2016).

The sizes of the molars and the width of the upper incisors (INCW) in Spalax sp. correspond to its small skull size; however, it has a significantly lower INCW/MSL ratio (14.8%, 13.6–16.1%; p = 0.013 to < 0.0001) than all other extant species (~17%, 13.7–21.4%; Table S5). In other words, this form is characterized by very narrow upper incisors. The width of the lower incisors appears to be slightly smaller than in S. microphthalmus (Table 3), whereas the longitudinal length of the incisors does not differ significantly.

Differences from S. microphthalmus (with low statistical significance) were detected in the lengths of m1 and m3, the width of m2, and the relative width of m3, which was slightly greater in Spalax sp. (Table 3). It should be noted that blind mole rat molars are rooted and that their crowns are strongly worn with age, resulting in high individual variability not only in size but also in the structure of the occlusal surface (Topachevsky 1969).

In Spalax sp., short glenoid cavities, combined with a relatively small difference in distance between their anterior and posterior edges (FGMIW and FGMAW, respectively), result in a rather sharp angle between the imaginary lines drawn along the surfaces of the glenoid cavities (Fig. S6A). A rough estimate of the angle between the glenoid cavities using only the FGMIW, FGMAW, and FGL measurements in Spalax sp. was approximately 18.8°, while in other species, this value exceeded 19° (19.7–23.8°).

In Spalax sp., the position of the upper alveolar ridge, which is not strongly developed, relative to the first upper molar differs markedly from that of S. microphthalmus (the alveolar ridge is located close to M1), but it is similar to that of the other members of the genus (Fig. 10).

A characteristic feature of Spalax sp. is the convex shape of the posterior margin of the hard palate. Moreover, this feature is detected already in young specimens with a non-overgrown suture between the basisphenoid and the basilar part of the occipital bone (Fig. 10). A distinct but small styloid process was observed in two specimens (S-218, S-224, TI), which are characteristic of species of the genus Nannospalax. The shape of the posterior palatine notch demonstrates individual variability. However, in Spalax sp., we did not find variants of the posterior palatine notch characteristic of S. microphthalmus (“straight edge”) or variants with a small notch in the anterior direction, as in the specimen of S. zemni (Fig. 10) or in S. graecus (a cleft in the area of the suture between the palatine bones: a trace of the styloid process (Topachevsky 1969)). Some specimens of S. giganteus also showed weak traces of the styloid process, but this feature is highly variable in this species.

Compared with S. microphthalmus, a notch in the sutura frontonasalis is usually present, but it is much smaller than in S. graecus (Fig. 11). The posterior edges of the nasal bones are usually slightly pointed rather than blunted, as in S. microphthalmus.

Age-related variability of the skull shape and progressive development of the sagittal and lambdoidal crests are typical of the blind mole rats, not only of Spalax (Topachevsky 1969). The lambdoidal crest of Spalax sp. is less powerful in comparison with this structure in S. microphthalmus or S. graecus. Skull height (CPSH) indirectly reflects the degree of development of this crest. The relative skull height in Spalax sp. is about 41.9% of skull length (38.0–43.6%). This value is comparable to the relative skull height observed in S. giganteus and S. uralensis, whereas the other species exhibit skull height indices exceeding 43% (Table S5).

The relative height of horizontal branch of the mandible of Spalax sp. was significantly (p = 0.0002) lower (MDH/MNDL: 27.9%, 25.9–31.2%) than in S. microphthalmus (29.7%, 23.7–33.5%; Table S5). The alveolar process without a ridge on the back surface that presents in S. graecus (Fig. 11B: 6) is much higher than the articular process in adult Spalax sp. The shape of the corono-alveolar notch (incisura corono-alveolaris) is analogous to that in S. graecus (Fig. 11B). The irregularity of the notch edge is due to the strong development of the anterior ridge of the alveolar process, which “comes into contact” with the coronal process (Topachevsky 1969). Thus, Spalax sp. represents the fourth known form with this type of corono-alveolar notch structure, together with S. graecus, S. arenarius, and the extinct S. minor (Topachevsky 1969).

The pattern of the molar crown masticatory surface in Spalax is very diverse on the one hand and not very different between species on the other. Among all forms, age-related variability predominates (Topachevsky 1969; Puzachenko 1991; see our Figs 3, 12). The most complex crown structure is observed in juveniles and sub-adult animals, up to 1–1.5 years of age. In these age groups, the occlusal surface of the molars clearly exhibits features typical of fossil Spalax species, including separation of the anterocone and paracone on M1; separation of the proto- and hypocone on all upper molars and of the proto- and hypoconid on the lower molars; and isolation of the protoconid on m1 (Topachevsky 1969; Sarica and Sen 2003). However, some discrepancy was observed between cranial age-related features (development of the sagittal and lambdoid crests and the shape of the parietal bones) and the degree of crown wear. In adult Spalax sp., as indicated by cranial characteristics, the tooth crowns were less worn than in adult S. microphthalmus or S. graecus specimens. In conclusion, Table 4 summarizes the main differences between the Central Caucasian Spalax sp. and the greater blind mole rat on skull measurements (see Table S5 for more details).

The karyotype of North Caucasian blind mole rats (2n = 62) does not differ significantly from other 62-chromosome karyotypes of Spalax species. All of these chromosome sets exhibit similar ratios of metacentrics, submetacentrics, and subtelocentrics with very similar morphology (Arslan et al. 2016). The sex chromosome morphology of Spalax species is also similar, with the exception of S. giganteus, in which the X chromosome is unequal-armed (submetacentric or subtelocentric).

S. microphthalmus (2n = 60) is the only species with a distinct chromosome set in the karyotypically conservative genus Spalax. We found no evidence of variation in chromosome number within this species. Throughout the sampled range of S. microphthalmus, including geographically distant populations from its eastern, western, northern, and southern limits, a consistent diploid number of 2n = 60 was recorded (L’apunova et al. 1974; Martynova 1977; Puzachenko and Baklushinskaya 1997; Arslan et al. 2016; this study). Compared with the 62-chromosome karyotype, the 60-chromosome set lacks one pair of chromosomes. We assumed that these are small submetacentric chromosomes, possibly homologous to the smallest pair of submetacentrics in the Spalax sp. set, based on the available G-banding data (Figs 4B, 5B). It is difficult to determine the location of the translocation of these elements in the 60-chromosome set. This requires the use of high-resolution staining methods. It is quite probable that the 62-chromosome North Caucasian blind mole rat has maintained the karyotype of the ancestral form common to itself and S. microphthalmus (2n = 60).

Analysis of SCs in pachytene spermatocytes of the 60- and 62-chromosome blind mole rats confirmed the existing karyological data. Notably, we observed an unusual “lumpy” configuration of the XY bivalent in both individuals, indicating that these karyotypes share this distinctive sex bivalent morphology. Because SCs have not been characterized in other Spalax species, direct intrageneric comparisons are not yet feasible. Nevertheless, the 60- and 62-chromosome blind mole rats are clearly similar in this respect. In contrast, species of the closely related genus Nannospalax exhibit a different sex bivalent configuration (Wahrman et al. 1985; Matveevsky et al. 2015, 2018, 2020). Therefore, a more comprehensive understanding of sex chromosome dynamics during meiosis in this group will require both additional material and broader taxonomic sampling.

The distribution area of Spalax sp. (2n = 62) and its boundary with S. microphthalmus (2n = 60) are not clearly defined. Karyotyped individuals with 2n = 62 were found in two localities in Kabardino-Balkaria (Dzuev and Shogenov 2004) and in the south of Stavropol Krai (our data) in marginal southern populations of S. microphthalmus sensu lato. The nearest localities where blind mole rats with 2n = 60 were found are the right bank of the Manych River and the left bank of the mouth of the Don River (our data), as well as the Taman Peninsula (Martynova 1977) located at a distance of 330–450 km from the localities of Spalax sp. (2n = 62). Whether a contact zone exists between the 60- and 62-chromosome mole rats, or whether they are allopatric, remains unresolved. It is worth noting that non-karyotyped animals morphologically similar to Spalax sp. (2n = 62) and different from S. microphthalmus (our data) were previously captured in the territory of Kabardino-Balkaria from the foothills to an altitude of about 2000 m a.s.l. (Dzuev and Shogenov 2004; Tembotova 2015).

Based on its position in the phylogenetic trees (Figs 7, S1) and the MJ haplotype network (Fig. 8), as well as high genetic distances (Table 2), Spalax sp. is clearly differentiated from all other Spalax species, including S. microphthalmus and S. giganteus, with which it was previously grouped (Dzuev and Shogenov 2003; Korobchenko and Zagorodniuk 2009). The cyt b sequence divergence between Spalax sp. and other species ranges from 8% to 12%, which substantially exceeds the 5% guideline typically applied for species-level differentiation in mammals (Baker and Bradley 2006). Our comprehensive genetic analysis reveals that Spalax sp. (2n = 62) represents a distinct phyletic lineage of blind mole rats with an independent evolutionary history, having undergone complete speciation and genetic differentiation from all other studied Spalax species.

The topology of the SpalaxmtDNA trees reveals an initial divergence of large-bodied blind mole rats into two distinct phyletic clades. The first clade comprises the western species, S. antiquus and S. graecus, whereas the second encompasses all remaining species distributed to the east. High statistical support for species branches combined with marginal support for deep nodes in ML trees (Fig. S1) indicate a rapid evolutionary radiation (explosive speciation) (Rokas et al. 2005) of Spalax, including S. microphthalmus and Spalax sp., from a common ancestral form (or some closely related forms) that occurred within a short geological time frame. The MJ haplotype network topology and the high number of mutational steps between species-specific haplogroups (Fig. 8) indicate prolonged independent evolution of Spalax species.

The divergence time of Spalax sp. and S. microphthalmus can be estimated by comparing our results with molecular dating available in the literature. Currently, two dated phylogenies of blind mole rats are available. The first is based on mitochondrial data and places the divergence between Spalax and Nannospalax at ~7.6 Mya (Hadid et al. 2012). The second is reconstructed using genome-wide SNP data and estimates this divergence at ~3.48 Mya (Yanchukov et al. 2026). Both studies focus on the phylogeny of Nannospalax and lack sufficient paleontological calibration points. Nevertheless, a comparative analysis of the tree topologies obtained in our study and those presented in the aforementioned publications allows us to extrapolate these dates to our results. Accordingly, the divergence between Spalax sp. and S. microphthalmus can be placed to the Early Pleistocene (Calabrian) based on the calibrations by Hadid et al. (2012), or to the Middle Pleistocene according to Yanchukov et al. (2026). Considering the geographic distribution of S. microphthalmus sensu lato, the isolation of these two forms appears to have been driven by one or more Ponto-Caspian transgressions, which periodically separated the Caucasus from the Eastern European Plain. In the first scenario, this event may correspond to the Apsheron Ponto-Caspian transgression, which lasted approximately 500,000 years (1.8–1.3 Ma; Forte and Cowgill 2013; Svitoch 2016). Alternatively, a Middle Pleistocene divergence between Spalax sp. and S. microphthalmus could have occurred during the Bakunian transgression, which lasted about 270 kyr (MIS 18–13; Krijgsman et al. 2019). The latter option aligns better with the morphological and paleontological evidence (see below).

The lack of an explicit phylogenetic structure of S. microphthalmus (2n = 60) may indicate the low efficiency of potential ecological and geographical barriers in the territory of its distribution. An alternative explanation is a recent rapid colonization of the area by a local ancestral population. The combination of high haplotype diversity (Hd) and low nucleotide diversity (π) suggests that the recent population may have originated from an ancestral population with a low effective population size. A recent population growth and range expansion is also indicated by the relatively high expansion coefficient (S/k) and the results of the neutrality tests, which showed negative values that were significant for D and Fs, but not significant for F*. The above data, combined with a low level of intraspecific genetic variability and the absence of geographic subdivision, suggest that the modern range of S. microphthalmus was formed as a result of rapid dispersal from one or several small refugia in historically recent times. Projecting this finding onto paleoclimatic events, we can hypothesize that during the Valdai glacial epoch, when forest and forest-steppe communities expanded across the present-day East European steppe and meadow biocenoses (Pakhomov 2006), the range of S. microphthalmus was confined to one or more refugia situated on the uplands of the East European Plain (such as the Central Russian or Volga uplands). The subsequent mid-Holocene warming resulted in climate aridization and the recovery of steppe landscapes, thereby providing the necessary conditions for the rapid expansion of the species throughout the Volga-Dnieper region. The abundance of Holocene fossils of S. microphthalmus (Topachevsky 1969) may serve as indirect evidence of this expansion.

The results of morphometric analysis and comparison of modern and fossil representatives of the genus Spalax allowed us to reveal a unique combination of qualitative and quantitative features in Spalax sp., distributed in the central part of the North Caucasus. In the skull structure of the Central Caucasian Spalax sp., some features stand out that may reflect an adaptation to digging, compared with S. microphthalmus, S. arenarius, etc.

According to our data (Table 3), the ratio of transverse to longitudinal length of the lower incisor is some higher in S. microphthalmus (105.1%) than in Spalax sp. (99.5%), but in this case we only consider this as a trend. According to (Topachevsky 1969), in S. giganteus, the index was 95.2–98.0 (mean)–108.1%. In the type series (Nogaisk = Prymorsk, Zaporizhia) of the Early-Middle Pleistocene S. minor this index was between 90.0–97.7–104.1% (Topachevsky 1969). The few specimens of S. minor from other Pleistocene localities from the lower Dnieper River and the Azov-Black Sea region had an index of 88.9–96.6% (Topachevsky 1969; Stadnik 2009). In the Middle Pleistocene form S. cf. microphthalmus, this index was within the range of 93.9–100.5–106.6% (Stadnik 2009), and, therefore, the latter form is closer to the recent S. microphthalmus than to Spalax sp. Thus, it is likely that Spalax sp. has the narrowest lower incisors among all the recent members of Spalax.

The average parameters of the molars in Spalax sp. differ significantly from those typical of S. minor (Stadnik 2009; see our Fig. 13). However, the minimum values for M1–M3 length and m1–m3 and the width of m2 and m3 fall within the range of S. minor. The variability in the length of M1 and M2, as well as in the length and width of M3, m1, m2, and m3 in the North Caucasian Spalax sp. covers the variability in the fossil form S. cf. microphthalmus, dated within a wide interval of ~0.76–0.37 Ma (Stadnik 2009). At the same time, the latter form exhibited, on average, narrower upper molars, a narrower m1, but a wider m2. For the lower molars, there is a “trend” of tooth size change from the Early Pleistocene S. minor to the modern S. microphthalmus. Molars from several Middle Pleistocene (~0.76–0.37 Ma) of S. cf. microphthalmus (Stadnik 2009), including the Otkaznoe section (0.42–0.37 Ma), and those of the recent Spalax sp. (2n = 62) occupy an intermediate position (Fig. 13).

A rough estimate of the angle between the left and right fossa glenoidea (Fig. S7A) in Spalax sp. from the Central Caucasus region was approximately 18.8°, a value lower than that observed in the other species. According to G.E Zubtsova (Zubtsova 1986), before blind mole rats bite the soil, the condyle of their condyloid process is within the middle third of the fossa glenoidea; and in this position of the mandible, the lower incisors are spread apart (Fig. S6B). The sharper angle between the glenoid cavities suggests a relatively small divergence of the upper parts of the lower incisors during digging, which reduces their overall efficiency (Zubtsova 1986). The biological interpretation could be as follows: a) Spalax sp. is more adapted to digging in dense soils; or b) this species generally has a relatively low adaptation to digging among all the other Spalax species.

The degree of development of the lambdoid crest is positively related to the development of the musculus rhomboideus capitis, which is attached along the entire length of this crest and terminates at the scapula (Gambaryan 1960). The strength of this muscle is related to the magnitude of the forces generated when the animal raises its head, which is of great importance in the digging process (ibid.). The development of the lambdoid crest is proportional to the load at the site of muscle attachment, and the crest itself increases with age, reaching its maximum size in adult and old males of all species. The width of the interorbital region and the development of the sagittal crest depend on the temporalis muscle, which is the primary retractor of the mandible. The function of this muscle is extremely important in the digging behavior of blind mole rats (Zubtsova 1986). During ontogeny, interorbital width and sagittal crest reach their minimum and maximum, respectively, in adult and old animals. In Spalax sp., the relatively low lambdoid crest, combined with a relatively wide interorbital region and weakly developed sagittal crest (Fig. 11), indicates a lower degree of specialization of the skull for digging activities compared with other Spalax species. It should be noted that the development of the crests appears to proceed more slowly during postnatal ontogeny in Spalax sp. than in S. microphthalmus. For this reason, sexually mature specimens of blind mole rats caught in Kabardino-Balkaria were often previously classified as subadults of S. microphthalmus.

The relatively narrow incisors of Spalax sp., compared with all other modern species, may also indicate lower digging efficiency. However, these same features can be interpreted as an adaptation to dense soils containing stone fragments or debris, which is typical of the known range of Spalax sp. in the foothill regions. It cannot be ruled out that the relatively high posterior part of the rostral region (increasing strength), compared with most other Spalax species, supports this latter interpretation.

We found that Spalax sp. has, on average, the smallest skull size among the extant species of the genus. If an increase in overall body and skull size is considered an evolutionary trend beginning with S. minor, the small size of Spalax sp. may be cautiously interpreted as an archaic character.

We also detected deviations from allometric relationships in many skull measurements, not only in comparison with S. microphthalmus but also relative to other species of the genus. Thus, the skull of these blind mole rats is not a simple “miniaturized copy” of the skull of the greater blind mole rat or of any other species, but instead exhibits a set of specific and original characters. To this should be added a more prolonged temporal development of skull structures directly associated with digging activity.

Summarizing the overall comparison of the Central Caucasian 62-chromosome form with other extant species of the genus Spalax, extinct S. minor and S. cf. microphthalmus, it can be concluded that, based on the examined complex of structural traits and certain quantitative parameters its skull and teeth retains several ancestral characters of the genus Spalax. Formally, this indicates a lower degree of adaptation to a subterranean lifestyle. Clarifying how these morphological features are reflected in the biology of this form will require filling a substantial gap in our knowledge of its ecology.

The probable origin of the described Spalax sp. can be traced back at least to the middle part of the Middle Pleistocene in Eastern Europe or Asia Minor. Considering the geographical position of the Otkaznoe locality, we assumed that the Spalax sp. found there is the direct ancestor of the modern Spalax sp. in the region. The available data do not exclude the possibility of a wider distribution of this form in the Middle Pleistocene. In our opinion, a revision of the fossil remains of “S. cf. microphthalmus” from the Middle–Late Pleistocene deposits is necessary.

It should be emphasized that the sample of Spalax sp. currently available in museum collections is relatively small and probably does not yet provide a complete picture of the extent of individual variability in particular skull measurements. For the same reasons, age-related variation and sexual dimorphism remain insufficiently studied.

Nevertheless, based on a combination of original morphological features and the available genetic data, we consider the Central Caucasian blind mole rat to represent a distinct species within the genus Spalax.

The combined data presented above suggest that the North Caucasian blind mole rat with karyotype 2n = 62, previously classified as the greater blind mole rat, S. microphthalmus (Dzuev 1989; Dzuev and Shogenov 2003; Dzuev et al. 2019) or S. cf. giganteus (Korobchenko and Zagorodniuk 2009), is a new species of the genus Spalax.

Class: Mammalia Linnaeus, 1758

Order: Rodentia Bowdich, 1821

Suborder: Myomorpha Brandt, 1855

Superfamily: Muroidea Illiger, 1811

Family: Spalacidae Gray, 1821

Subfamily: Spalacinae Gray, 1821

Genus: Spalax Güldenstädt, 1770