Erebia manto is a species of the sub-alpine zone, predominantly present around or just below the tree line. It is associated with humid habitats with tall growing grasses, both natural and cultivated (Warren 1936; Sonderegger 2005). Its distribution is restricted to Central Europe, where it inhabits the three major massifs: the Alps, the Pyrenees and the Carpathians, and some smaller ones: the Massif Central, the Vosges and the Jura (Warren 1936; Willien 1985; Lafranchis 2014; Székely 2008). In the Balkan Peninsula some isolated populations exist in Bosnia (Jaksic 1988; Varga 2014). Records of the Durmitor in Montenegro (van Swaay et al. 2010) are erroneous (Jaksic 1988; Dug 2013; Varga 2014). It is absent from the Italian and the Iberian Peninsula south of the Pyrenees. It has been reported from the Picos de Europa (Gómez Bustillo and Fernández Rubio 1974; Mortera 1996), but the reliability of this observation has been doubted (González Fernández 2010; Mortera 2013). Erebia manto was originally described from the eastern part of the Alps, but Denis and Schiffermüller (1775) fixed no holotype, nor a type locality. Delimitation of the nominotypical subspecies from ssp. mantoides, the latter of which is widespread in the Alps, is not clear, neither morphologically nor geographically. For practical reasons, we refer to all our Alpine samples as ssp. mantoides.

Figures 1, 2. 

Morphotypes and genetic lineages described for E. manto. 1. Valve shape of typical representatives of the three morphotypes after Cupedo (1997). A manto type, Col de la Croix (CH-VD), prep.nr. 1503. B bubastis type, Col de la Madeleine (F-73), prep.nr. 1113. C vogesiaca type, Rainkopf (F-68), prep.nr. 1705. Scale bar: 1 mm. 2. Neighbour-joining phenogram by Schmitt et al. (2014), based on Nei’s genetic distances (Nei 1972). Numbers represent individual samples. Bootstrap values >50% (1000 iterations) are given at the nodes.

A total of 152 E. manto adult individuals were collected at 15 localities between 2003 and 2013 (Table 1 and Fig. 3). Because of their geographic isolation, we treat the different localities as populations. As an outgroup, five males were sampled of the closely related E. eriphyle (Freyer, 1836). Netted individuals were paralysed with ethyl acetate vapour, and two legs were transferred to 96% ethanol before death. Specimens included in this study are all registered on the Barcoding of Life Database (Ratnasingham and Hebert 2007) with detailed information on collecting localities and identification (dataset DS-manto: https://doi.org/10.5883/DS-MANTO). Vouchers are deposited in the collection of Naturalis Biodiversity Center (RMNH; Leiden, the Netherlands).

Figure 3. 

Sampling locations in the Alps, the Jura and the Vosges. Light grey: mountains > 1000 m. Dark grey: lakes and coastline. Sample numbers refer to Table 1.

Genomic DNA was extracted from the legs, or in a few cases from part of the abdomen, with a Macherey-Nagel NucleoMag 96 Tissue magnetic bead kit on a Thermo Fisher KingFisher flex system. Polymerase Chain Reaction (PCR) was used to amplify a fragment of 1,496 base-pairs (bp) of the cytochrome oxidase subunit I gene (COI) from the mitochondrial genome, in two parts. The first part of 658 bp, near the 5’-P side of the gene, also known as the DNA barcode region (Hebert et al. 2003) was amplified using a 1:1 primer-mix of LepF1 (5’ ATT CAA CCA ATC ATA AAG ATA TTG G 3’) and LCO1490 (5’ GGT CAA CAA ATC ATA AAG ATA TTG G 3’) as forward primers, and a 1:1 primer-mix of LepR1 (5’ TAA ACT TCT GGA TGT CCA AAA AAT CA 3’) and HCO2198 (5’ TAA ACT TCA GGG TGA CCA AAA AAT CA 3’) as reverse primers (Folmer et al. 1994; Doorenweerd et al. 2014). The second part that we amplified, from here on indicated as the 3’-P fragment, has a 32 bp overlap with the 5’-P fragment, and a target length of 870 bp. This fragment was amplified with forward primer SeqIntCOIf (5’CWT CWT TTT TTG AYC CAG CWG GAG 3’) and reverse primer LepLEUr (5’ CCA TTA CWT ATA RTC TGC CAT ATT 3’) (Vila and Björklund 2004). M13 forward (5’ TGT AAA ACG ACG GCC AGT 3’) or reverse (5’ CAG GAA ACA GCT ATG AC 3’) tails were attached to the forward and reverse primers respectively, and used for sequencing. Each PCR reaction included 40 cycles, consisting of 3 minutes initial denaturation at 94 °C, 15 seconds cycle denaturation at 94 °C, 30 seconds cycle annealing at 50 °C, and 40 seconds cycle extension at 72 °C. A final extension at 72 °C for 5 minutes completed the reactions. Bidirectional Sanger sequencing was outsourced to BaseClear, Leiden, The Netherlands. The resulting chromatograms were checked for quality and congruence in Geneious R6.1.8 and the resulting sequences were managed using VoSeq 1.7.4 (Peña and Malm 2012) prior to uploading to BOLD and Genbank (accessions MH670083–MH670241). Sequences were trimmed to remove any missing data, and specimens with too much missing data were removed; leaving 133 E. manto samples for haplotype analyses. A median-joining haplotype network was constructed using Network v4.6.0.0 (Fluxus-Engineering) and optimized for publication in Adobe Illustrator CS6.

An estimate of the chronology of differentiation phases was obtained through StarBEAST v2.5 (Heled and Drummond 2010), which allows for coalescent histories of the clusters of specimens. We applied a HKY substitution model with four gamma categories, a constant population function, mitochondrial ploidy, a birth-death species tree and ran the analyses for 40 million generations. Sufficient sampling of all parameters was checked using Tracer 1.7. To calibrate our phylogeny to absolute ages, we applied an external calibration based on a four-gene Erebia-wide phylogeny (Peña et al. 2015). That study in itself also relied on an external calibration, with date-ranges derived from a phylogeny of Nymphalidae that was calibrated using fossils and host plant family ages (Wahlberg et al. 2009). Although this double-external calibration undoubtedly increased the margin of error, the 95% highest posterior density ranges of the node age estimates were used as normally distributed priors, and the results should roughly equate to the same estimate ranges as in the primary source. It is important to note that if future studies prove the estimates in the primary work to be inaccurate, this cascades down to the study of Peña et al. (2015) and the study presented here. Nonetheless, we felt that this was the most accurate dating method currently available; more reliable than using mutation rates for COI derived from other groups of insects. Resulting trees were compiled using TreeAnnotator with a 25% burn-in, visualized with FigTree and prepared for publication using PaintShop Pro X6.

A COI sequence of 1,496 bp was obtained for 133 individuals (Table 1). Twenty-nine haplotypes were found, defined by forty-seven variable sites, accounting for 3.14% of the sequence. Of these, 46 had one substitution, at one site two mutations were detected. Thirty-four substitutions (71%) were transitions. Eight single-nucleotide polymorphisms (SNP’s) were singletons, 40 were parsimony informative, 18 of which were located on the 5’-P fragment and 22 on the 3’-P fragment. The percentages of parsimony informative SNP’s in the 5’-P and the 3’-P fragment were 2.74 and 2.63 respectively, in contrast with what we expected on the basis of Lunt et al (1986). The median-joining haplotype network (Fig. 4) is highly differentiating: 26 out of 29 haplotypes (90%) were private to a single population. The network shows five well differentiated clusters. Each morphotype is represented by one or two clusters, with the exception of two haplotypes that were shared by populations belonging to different morphologic groups (Table 2). Cluster A includes the haplotypes of the Pyrenean ssp. constans. In cluster B the samples of the ssp. mantoides, roberti and osmanica are found, and one sample of ssp. vogesiaca. Cluster C covers the ssp. bubastis and willieni, and contains one sample of ssp. mantoides. The last two clusters are confined to one subspecies each: ssp. vogesiaca (D) and ssp. trajanus (E). The cluster arrangement is star-like, with cluster B in the centre. Cluster A is 16 mutational steps apart (1.1% pairwise difference), and C, D and E three steps (0.2%). The outgroup Erebia eriphyle is connected to the Pyrenean population by 89 SNP’s (5.9% pairwise difference).

Figure 4. 

Median-joining network of the COI haplotypes obtained from 133 individuals of Erebia manto. Circle diameter is proportional to the frequency of the haplotype. Black dots represent missing haplotypes. Clusters are indicated by colours and are designated A through E. Deviant colours within a cluster indicate samples with deviant morphotype.

The estimates of divergence times of the haplotype clusters shown in Fig. 4, based on the StarBEAST results, are shown in Fig. 7. Following the estimates from the concatenated dataset analyses in Peña et al. (2015; fig. S8), we used a normal distribution of 9.2–14.5 (5%–95%) to date the split between E. manto and E. eriphyle, and 1.8–4.5 (5%–95%) to date the split between E. manto constans and all other groups. Although their per-gene partitioned analyses (Peña et al. 2015; fig. S7) estimated significantly younger ages for these splits (4.9–10.0 and 0.9–3.2, respectively), the authors suggest that this is due to difficulties of fitting a substitution model with genes that have saturation (i.e. COI) and that the concatenated dataset ages are likely more accurate. Applying these calibration points to our data places the intraspecific divergence of E. manto in the Pleistocene. Pyrenean E. manto constans splits off from all others most likely in the early Pleistocene. All other groups diverge from each other in the second half of the Pleistocene, but their relationships and relative timing of divergence remains statistically unsupported with our data.

In the absence of high-resolution genomic data, as is the case here, the most reliable reconstruction of phylogenetic relationships is achieved by an integrative ‘total evidence’ approach, including morphological data (Simonsen et al. 2006; Ferri et al. 2009; Pazhenkova and Lukhtanov 2018). Unlike adaptive traits like wing pattern (Brown et al. 2000; Simonsen et al. 2006), exoskeletal structures like wing venation and genitalic appendages have been proven phylogenetically relevant (Damgaard and Sperling 2001; Pisani et al. 2007; Schlick-Steiner 2010; Shi et al. 2015). We discuss the three datasets integratively, addressing incongruencies, and combine them to the greatest possible extent.

Cluster A. The deep genetic split between the Pyrenean population and the rest of the network is in line with the allozyme results (Schmitt et al. 2014). The allozyme study revealed that this lineage includes the populations from the French Massif Central (Schmitt et al. 2014). A cohesion of Pyrenean and Massif Central populations is not unique in the genus Erebia (Hinojosa et al. 2018). Morphologically, the clusters A and B share the manto type valva. Evidently, the Pyrenees harbour a well-defined genetic lineage, not recognised morphologically. Populations of cluster A exhibit a strongly melanistic wing coloration which, however, we do not consider reflecting the genetic uniqueness of this group (Schmitt et al. 2014), as melanism has a high adaptive value (Roland 1982; Guppy 1986; Clusella Trullas et al. 2007, 2008), and the coexistence of melanistic and richly decorated populations in the same (intraspecific) genetic lineage has been evidenced at least in E. epiphron (Schmitt et al. 2006) and E. euryale (Schmitt and Haubrich 2008). We will refer to cluster A as the constans lineage.

Cluster B covers the remaining manto type populations. Geographically, this lineage occupies the Alps (ssp. manto and ssp. mantoides), the Tatra and the Fatra (ssp. roberti) and a small area in Bosnia-Herzegovina (ssp. osmanica). We will refer to it as the manto lineage. In the Alpine populations of this cluster, the allozyme based neighbour-joining tree reveals three sub-clusters (Fig. 2). In the COI network, only the sub-cluster from the south-eastern Alps is separated. It is represented by the Dolomites, Hohe Tauern and Bosnia samples. The other two clusters are merely represented by a single sample each.

Cluster C in the network coincides with the bubastis group as a whole. This group, hitherto only characterised morphologically, evidently constitutes a separate genetic lineage, different from the remaining E. manto populations in this study, and hereafter referred to as the bubastis lineage. Its distribution area is strongly disjunct: it splits up into three, geographically widely separated sub-areas in Switzerland, France and Italy (Fig. 6). Morphologic differences are weak and phylogenetically irrelevant (Cupedo 1997), but their spatial distributions, at both sides of the Alpine main chain, indicate that the extant populations trace back to at least two glacial refugia. Hence their separation must be at least of Würm age (roughly 10–70 Ka BP). It is logically inferred that the separation of the bubastis lineage as a whole must be of at least Riss age (roughly 120–180 Ka BP). This cluster is the only one in which a third, even younger (i.e. post Würm) level of morphological differentiation has been described. Both in the French and in the Swiss sub-area, neighbouring, but nowadays isolated populations differ in the allele frequency for a conspicuous wing pattern element (Sonderegger 1983; Cupedo 1997), an apparent case of postglacial genetic drift. This observation supports the Würm glacial age of the disjunction of these sub-areas. In cluster C appears one population belonging to the manto morphotype: the Roselend sample (manto type) shares its haplotype (No 18) with the nearby Pralognan sample (bubastis type). Both the complex of morphological differences and the otherwise strict separation of both valve types, both in the network and geographically, make a close relationship improbable. Their geographic vicinity is an argument in favour of introgression. We hypothesise a historic unidirectional introgression into the Roselend population.

Cluster D contains the haplotypes of the Vosges sample of ssp. vogesiaca. Cluster E includes the haplotypes of ssp. trajanus. They are discussed in combination, because their independent origin appears not unquestionable when considering all available data. In the haplotype network these clusters constitute two independent branches, rooting in cluster B, and at a distance of three snp’s from it. In the allozyme based N-J tree they are seen as two branches of the same genetic lineage (Fig. 2). Schmitt et al. (2014) reported, though, that not all their analyses were concordant. The close relation suggested by the N-J tree was supported by the STRUCTURE analysis, but was not confirmed by their allele analysis. Evaluating all data, they concluded that the clustering in the N-J phenogram is an artefact, caused by genetic depletion. Besides, they consider the geographic distance arguing against derivation of a common ancestor. Morphologically, though, both subspecies share a valve type that is not found elsewhere, strongly suggesting a synapomorphy. Altogether, the available data are insufficiently consistent to establish unambiguously the relation between ssp. vogesiaca and ssp. trajanus. In order to avoid premature splitting, we provisionally consider them one entity, further on referred to as the vogesiaca lineage. The Jura population, belonging to this morphotype, shares its single haplotype with the Rofan population (manto type). Here too, a close relationship of the Jura and Rofan populations is improbable because of the pronounced genital differences. All populations sharing the Jura genital structure are found in haplotype clusters separated from cluster B by three mutational steps. Asymmetric introgression is regarded the most probable cause of this kind of discrepancies (Stephanović et al. 2016). Given the considerable geographic distance, a recent involvement of the actual Rofan population can be excluded. We suppose a historic introgression, probably far back in time.

Past disjunctions led to today’s biogeographic patterns, and vice versa: present distribution patterns yield keys to the reconstruction of past disjunctions or coherences. It is in this light that we assess three questions: 1) Is the coherence of the bubastis lineage, supported by mtDNA and morphology, compatible with its highly disjunct distribution pattern? 2) Do the highly disjunct geographic distributions of ssp. vogesiaca and ssp. trajanus argue against a genetic coherence, as suggested by morphology and the allozyme based NJ tree? And 3) is the supposed historic introgression of mtDNA into the Roselend population biogeographically plausible?

Ad 1. The populations of the bubastis type are concentrated in three, disjunct sub-areas in Switzerland (Sonderegger 1983, 2005), Italy (Floriani 1965; Baldizzone 1971) and France (Cupedo 1991, 1997) (Fig. 6). Together, they are restricted to the southern part of the western Alps, and sharply separated from populations of the manto type, which occupy the northern massifs. Both in France and in Switzerland the dividing line between the manto type and the bubastis type populations coincides with a well-known biogeographic boundary in the genus Erebia, separating species or subspecies of northern/north-western alpine and southern alpine origin. In the Swiss canton of Wallis it basically follows the Rhone valley east of Martigny. Upstream, however, some taxa originating from the southern alpine border crossed the Rhone, and dispersed westward on its right flank (Schweizerischer Bund für Naturschutz 1987). Taxa of northern or western origin colonised the northern Rhone bank from the west. Halfway they met, building a zone of encounter or admixture (Fig. 5, 6). In the genus Erebia, the northern group is represented by E. euryale isarica, E. epiphron aetheria, E. pronoe vergy, E. oeme and E. manto mantoides. In the southern group we find E. euryale adyte, E. epiphron valesiana, E. pronoe psathura, E. alberganus, E. mnestra and E. manto bubastis. For detailed distributional data and charts see Gonseth (1987), Schweizerischer Bund für Naturschutz (1987), and Sonderegger (2005). The distributions of the manto type and bubastis type populations fit exactly in this scheme. They come into contact just west of Leukerbad. The Tessin population testifies of the former connection to the southern Alpine foothill. In the French Alps the border between southern and northern taxa is equally sharp and follows the Isère valley. It separates the distribution areas of the same taxa as in Wallis (Willien 1985; Lafranchis 2014; Diatheo 2016), except E. epiphron, for which adequate morphologic data are missing in this region, and E. pronoe psathura, which is absent here. The northern taxa occupy the French Pre-Alps too, but are absent from the left bank of the Isère north of Grenoble (Fig. 5). Here too, some taxa of southern origin crossed the upper Isère and invaded the Beaufortin in the Mont Blanc massif (e.g. E. euryale adyte, see Fig. 6). The bubastis type and the manto type populations strictly respect this borderline. As in the Swiss part of their distribution, they are a textbook example of a sibling pair with a north-south disjunction. Altogether, the area of the bubastis morphotype as a whole fits zoogeographic dogmas. It covers the southern Alps from Lake Como in the east to the Provence in the west. This part of the Alps has been established as a glacial refuge and genetic hotspot for many plant and animal species (Schönswetter et al. 2005; Pauls et al. 2008; Schmitt and Haubrich 2008; Louy et al. 2014). The morphological and molecular coherence of the extant populations, combined with their discontinuous distribution pattern are highly suggestive of relict populations in a formerly contiguous distribution area. The manto cluster and the bubastis cluster are distinct genetic lineages, apparently resulting from an ancient (Riss or pre-Riss) north-south disjunction in the western Alps.

Figures 5, 6. 

Distribution areas of some Erebia taxa in the western Alps. Dashed line: outskirt of the Alps, black: rivers and coastline, grey: lakes and sea, RR – river Rhone, RI – river Isère, MS – Mediterranean sea, LG – Lake Geneva, LCo – Lake Constance, LM – Lake Maggiore, LC – Lake Como. 5. Species originating from northern alpine refugia. Red: E. manto mantoides, blue: Erebia oeme. 6. Species originating from southern alpine refugia. Red: E. mnestra, blue: E. euryale adyte. Red dots: E. manto, bubastis morphotype. Distribution data after Willien (1985), Cupedo (1997), Schweizerischer Bund für Naturschutz (1987), Sonderegger (2005), Lafranchis (2014) and Diatheo (2016).

Ad 2. The distribution pattern of the populations exhibiting the vogesiaca morphotype is extremely disjunct, with representatives in the Vosges and Jura, and in the southern Carpathians. This mountain chain borders on the Pannonian basin, an area that has been an important glacial refuge for a wide variety of organisms (Willis and van Andel 2004; Feurdan 2005; Janovska and Pokorný 2008; Kuneš et al. 2008). A westward expansion from this refuge, north of the Alps, is paralleled in a broad range of organisms, including butterflies. It is considered one of the main paradigms of postglacial recolonization (Hewitt 1999, 2000; Schmitt 2007). It has been postulated that Erebia euryale, a species with the same altitudinal preference as E. manto, followed this post-glacial migration route as well (Schmitt and Haubrich 2008). In this light, the presence of genetically related populations in the southern Carpathians and in the Vosges and Jura is not surprising. The presence of manto type populations in between (in the Fatra and the Tatra) does not argue against this scenario, because the dispersion of the vogesiaca type populations preceded the colonization of Tatra and Fatra. Schmitt et al. (2014) presume at least a Riss age of the separation of ssp. vogesiaca, and a post Würm age of the split between Slovakian and northern Alpine populations. We conclude that derivation from a common ancestor of ssp. vogesiaca and ssp. trajanus is not geographically implausible.

Ad 3. We suppose a historic introgression of a bubastis type population (sample 11) into a manto type population (sample 4) in the Beaufortin massif (Mont Blanc group, France). In this region, the Isère valley marks the border between both types. They meet upstream, somewhere between sampling sites 4 and 11. This makes a historic encounter and intermingling plausible. Nowadays, E. manto occurs in localised populations, and the individual range of action hardly exceeds the population boundaries. Genetic exchange most probably goes back to the recolonization phase in the early post-glacial, when bio zones were lower and hence more extended, and distribution ranges were shifting continuously.

As to the timing of the nodes in the graph, our data are merely a preliminary assessment, due to the large confidence intervals and the lack of support after the third node, but they do set some boundaries within which differentiation likely occurred. Subspecies constans was the first lineage to diverge from all remaining E. manto clusters. Its branch is well supported, and the split dates from the Pliocene/Pleistocene boundary. The younger nodes are probably associated with the Pleistocene glaciations. Branches are poorly supported, though, and no conclusions can be drawn on the sequence of branching. At this point the allozyme results showed more resolving power. They support two differentiation phases in this part of the time scale. The differentiation of the manto lineage, the ssp. vogesiaca and the ssp. trajanus was considered to be associated with Riss (MIS6) or even older. The same is argued in this paper for the bubastis lineage, judging from its nested morphologic structuring and its geographic distribution pattern. The split of the three Alpine subgroups of the manto lineage, and the differentiation of ssp. roberti in the Tatra, were proven considerably younger, presumably of Würm age (Schmitt et al. 2014). In summary, the phylogenetic history of E. manto seems to have gone through three main stages: the separation of the constans lineage around the Pliocene/Pleistocene boundary, the radiation of the manto and bubastis lineages and the ssp. vogesiaca and trajanus during the earlier glaciations, and differentiations within these lineages during Würm.

Figure 7. 

Bayesian tree of the clusters A-E in Fig. 4 and E. eriphyle as outgroup, constructed with StarBEAST2. Support values in brackets; blue vertical bars indicate 95% highest posterior density range; dotted red lines indicates the Pliocene-Pleistocene boundary and the Miocene-Pliocene boundary, after Gibbard and van Kolfschoten (2005) and Deutsche Stratigraphische Kommission (2016); grey zone marks the period of the Pleistocene glaciations, i.e. MIS2 – MIS16 according to Ehlers and Gibbard (2007).

Intraspecific differentiations, leading to genetically and morphologically different but geographically coherent groups of populations, are generally considered the result of incomplete vicariant speciation processes (Mayr 1963; Coyne and Orr 2004). Because, in the case of E. manto, these groups live in complete allopatry, the degree of their reproductive compatibility remains unknown. It has repeatedly been suggested that the E. manto complex might actually comprise two sibling species (Cupedo 1997; Schmitt et al. 2014). Morphologically, the differences in the male genitalia between the manto type and the bubastis type are such as, in the genus Erebia, only found between validated species (Fig. 1), suggesting a specific status of the bubastis lineage (Cupedo 1997). In the allozyme study, the genetic distance between the Vosges population and the remaining populations was almost as large as between E. manto and E. eriphyle, suggesting a specific status of the ssp. vogesiaca (Schmitt et al. 2014). In both cases it was emphasized that definite proof is lacking. Indeed, in the two cases of introgression postulated above, exactly the taxa under discussion were involved. That implies that the three groups were not, at least not completely, reproductively isolated from each other at the moment of introgression, a situation that most probably persists today. Consequently, their taxonomic status should be considered intraspecific. As to the constans lineage, its complete allopatry prevents any observation on the degree of reproductive compatibility to other groups. Altogether, any attempt to split up E. manto into two species (Ziegler 2005–2016; Centre Suisse de Cartographie de la Faune 2016; Lepiforum 2016) is at least premature. Nonetheless, considering all described subspecies equivalent does not do justice to their nested structure and the different ages of their origin. Nomenclature can only reflect a hierarchic phylogenetic structure of a species by pooling subspecies into groups, according to art. 6.2 of the code (ICZN 2000). Four such groups can be recognised: Erebia manto (group manto), encompassing the ssp. manto, mantoides, praeclara and osmanica; E. manto (group constans), consisting of ssp. constans; E. manto (group bubastis), containing the ssp. bubastis, valmaritima and willieni, and E. manto (group vogesiaca), comprising at least ssp. vogesiaca. The southern Carpathian ssp. trajanus is provisionally placed in this group, for lack of consensus between the available data regarding its relation to ssp. vogesiaca.