Extensive variation in chromosome number and genome size in sexual and parthenogenetic species of the jumping-bristletail genus Machilis (Archaeognatha)

Parthenogenesis in animals is often associated with polyploidy and restriction to extreme habitats or recently deglaciated areas. It has been hypothesized that benefits conferred by asexual reproduction and polyploidy are essential for colonizing these habitats. However, while evolutionary routes to parthenogenesis are manifold, study systems including polyploids are scarce in arthropods. The jumping-bristletail genus Machilis (Insecta: Archaeognatha) includes both sexual and parthenogenetic species, and recently, the occurrence of polyploidy has been postulated. Here, we applied flow cytometry, karyotyping, and mitochondrial DNA sequencing to three sexual and five putatively parthenogenetic Eastern-Alpine Machilis species to investigate whether (1) parthenogenesis originated once or multiply and (2) whether parthenogenesis is strictly associated with polyploidy. The mitochondrial phylogeny revealed that parthenogenesis evolved at least five times independently among Eastern-Alpine representatives of this genus. One parthenogenetic species was exclusively triploid, while a second consisted of both diploid and triploid populations. The three other parthenogenetic species and all sexual species were diploid. Our results thus indicate that polyploidy can co-occur with parthenogenesis, but that it was not mandatory for the emergence of parthenogenesis in Machilis. Overall, we found a weak negative correlation of monoploid genome size (Cx) and chromosome base number (x), and this connection is stronger among parthenogenetic species alone. Likewise, monoploid genome size decreased with elevation, and we therefore hypothesize that genome downsizing could have been crucial for the persistence of alpine Machilis species. Finally, we discuss the evolutionary consequences of intraspecific chromosomal rearrangements and the presence of B chromosomes. In doing so, we highlight the potential of Alpine Machilis species for research on chromosomal and genome-size alterations during speciation.


Introduction
In animals, sexual reproduction prevails and asexuality is uncommon in most taxonomic groups (Suomaleinen et al. 1987). It is assumed that asexual organisms have originated from sexual ancestors (Bell 1982), and this transition can occur multiple times independently within a taxonomic group (Stenberg et al. 2003;Schwander and Crespi 2009;Bode et al. 2010;Elzinga et al. 2013). Apart from some entirely asexual clades (e.g., bdelloid rotifers or darwinulid ostracods, the so-called ancient asexuals; Judson and Normark 1996), strictly asexual animal species are rare (Bengtsson 2009;Vrijenhoek and Parker 2009). Instead, asexuality often occurs as alternative reproductive strategy within a sexual taxon. Vandel (1928) first used the term geographic parthenogenesis to describe the phenomenon of distinct geographic distribution of sexual and asexual populations within a species. Since then, many examples have been reported of parthenogenetic lineages being more common at higher altitudes, higher latitudes, and in anthropogenically disturbed or extreme habitats (e.g., recently deglaciated areas), compared with their sexual relatives. Following this, the advantage of asexuals over sexual congeners for colonizing new habitats has been attributed to their doubled reproductive potential (as predicted by the twofold cost of sex (Maynard Smith 1978)) and to highly specialized ("frozen niche variation model") or generalized ("general purpose model") genotypes (Vrijenhoek and Parker 2009). It has been unclear, although, to what extent asexuality itself (Cuellar 1977;Glesener and Tilman 1978;Law and Crespi 2002;Maniatsi et al. 2011) or frequent correlates like polyploidy (Zhang and Lefcort 1991;Stenberg et al. 2003;Comai 2005;Adolfsson et al. 2010) or hybridization (Kearney 2003(Kearney , 2005Ghiselli et al. 2007) promote the spread of organisms into these habitats. The proposed advantages of polyploidy and hybridization include protection against deleterious mutations (Comai 2005) and increased heterozygosity (Kearney 2005), respectively. However, in the latest review of empirical data (Lundmark and Saura 2006), polyploidy turned out to be the factor best explaining distributional patterns in species containing asexual populations.
Parthenogenesis is mostly associated with polyploidy in plants, and, to a lesser extent, this correlation holds for animals as well (Suomalainen 1940;Otto and Whitton 2000;Choleva and Janko 2013). The low occurrence of polyploid lineages in animals has been explained by dioecy (i.e., male and female gametes come from different individuals) and chromosomal sex determination (for a detailed discussion, see Otto and Whitton 2000 and references therein), but in fact, polyploidy may have been simply overlooked in many animal species (Mable 2004). Polyploidy in animals is thought to emerge more likely in parthenogenetic diploid populations, while in plants, polyploidy usually predates asexuality (Otto and Whitton 2000). Already Bell (1982) suggested that in animals, polyploidy may be fundamental for the persistence of parthenogenetic lineages but not essential for their emergence. Despite these theoretical predictions, empirical studies on the interactions between asexuality and polyploidy and their ultimate effect on the evolution and spatial distribution of species remain scarce.
To enhance our understanding of the interrelation between asexuality and polyploidy, it is thus fundamental to collect empirical data from previously neglected animal groups with high incidence of parthenogenesis. Jumping bristletails (Insecta: Archaeognatha) of the genus Machilis meet these criteria. Within the Eastern Alps, 25 nominal species are known and at least nine of them putatively reproduce via parthenogenesis, as only females have been reported (Wygodzinsky 1941;Janetschek 1954). Moreover, Dejaco et al. (2012) hypothesized the occurrence of polyploidy in a study including three Machilis species.
In this study, we apply relative genome-size measurements, karyotyping, and mitochondrial DNA sequencing to sexual and putative parthenogenetic Machilis species to address the following questions: (1) Did parthenogenesis arise once or multiple times independently across species?

Specimen collection
We focused on eight species of which three reproduce sexually (M. helleri, M. hrabei, and M. lehnhoferi) and five putatively via parthenogenesis (M. fuscistylis, M. pallida, M. engiadina, M. ticinensis, and M. tirolensis). In the parthenogenetic species M. ticinensis, males have been reported from the southern Swiss Alps (Wygodzinsky 1941), potentially indicating a pattern of geographic parthenogenesis. In the absence of species-specific multilocus markers to prove asexual reproduction, we base our definition of parthenogenetic species on the following observations: (1) Four species, M. engiadina, M. fuscistylis, M. pallida, and M. tirolensis, have been collected numerously over the past 60 years but males were never found.
(2) Between 2009 and 2013, we collected over 2000 specimens from these and other parthenogenetic Machilis species throughout the Eastern Alps. In doing so, we considerably widened the known distribution ranges of these four species (e.g., Rinnhofer et al. 2012;T. Dejaco, unpubl. data) but again did not find males. (3) In the sexual species we sampled, males occur at roughly even ratios and males do not differ from females in their habitat preferences. It is thus unlikely that we overlooked males in presumably parthenogenetic species. (4) The fifth parthenogenetic species included in this study, M. ticinensis, was originally described as sexual species from the southern Swiss Alps (Wygodzinsky 1941). Between 2009 and 2013, we repeatedly sampled female M. ticinensis on the northern side of the Alps (Austria and Switzerland), but never encountered male specimens. In contrast, male and female M. ticinensis were found at an even ratio in one sampling location in the southern Swiss Alps (MIR). In total, 209 specimens were sampled at 41 geographically representative sites throughout the known distribution range of each species (Fig. 1; online supplementary material Table S1). Moreover, 12 specimens of four additional Machilis species were sampled: M. glacialis, M. inermis, M. mesolcinensis, and one unknown species, which morphologically keyed out as M. glacialis, but was named M. sp. A due to its distant position from M. glacialis in the mitochondrial phylogeny. These additional species were included in our phylogeny but excluded from other analyses due to poor geographic coverage. One specimen of Lepismachilis y-signata was used as out-group, resulting in a total of 222 specimens. Species were determined using the identification key in Palissa (1964) and original species descriptions (Riezler 1941;Wygodzinsky 1941;Kratochv ıl 1945;Janetschek 1949). Sexing of individuals was based on secondary sexual organsfemales are easily recognized by the presence of an ovipositor, detectable via screening a specimen's external morphology. In some cases, specimens were kept alive in a climate cabinet at 19°C until further use.

DNA extraction, PCR conditions, and phylogenetic reconstruction
Genomic DNA was extracted from muscle tissue using the Mammalian Genomic DNA Miniprep Kit (Sigma-Aldrich, St. Louis, MO). Two newly designed primers, MachF5 (5 0 -TAGTTATACCYATYATAATYGGHGG-3 0 ) and MachR7 (5 0 -CCTATRATAGCAAATACTGCYCC-3 0 ), were used to amplify approx. 750 bp of the mitochondrial cytochrome c oxidase 1 gene (CO1). PCR conditions were as follows: 95°C for 2 min, 35 cycles (94°C for 30 sec, 50°C for 45 sec, 72°C for 90 sec), and 72°C for 10 min. Amplicons were checked via gel electrophoresis and purified using a mastermix containing the enzymes Exo1 (1 U/lL) and FastAP (0.05 U/lL) (Thermo Fisher Scientific Inc., Waltham, MA) and applying an incubation step at 37°C for 15 min, followed by denaturation of enzymes at 80°C for 15 min. Sanger sequencing was conducted by a commercial sequencing facility (Eurofins MWG Operon, Munich, Germany) using the forward primer MachF5. All sequences were deposited in GenBank (accession numbers KJ501697 -KJ501918). Sequences were checked for correct amino acid translation and aligned using the Clu-stalW algorithm implemented in MEGA 5 (Tamura et al. 2011), yielding a final alignment of 703 bp. The HKY+G model of nucleotide substitution was determined as best fitting our data based on the Akaike information criterion using jModeltest 2 (Guindon and Gascuel 2003;Darriba et al. 2012). Net between-and within-species genetic p-distances were calculated in MEGA 5. The Mann-Whitney U-test was used to test whether within-species p-distances differed significantly between sexual and parthenogenetic species.
A Bayesian inference tree was constructed with MrBayes 3.2 (Ronquist et al. 2012). Three partitions were specified corresponding to codon positions. Two parallel runs, each consisting of three heated and one cold chain, were run for 10 6 generations, sampling trees every 1000 generations. Convergence was checked using the average standard deviation of splits frequencies, which fell below 0.01 after approx. 8.5 9 10 5 generations. The first 500 trees were discarded as burn-in. White and black symbols correspond to sexual and parthenogenetic species, respectively. Whenever more than one species was found at a location, species-specific symbols are connected by lines with that location.

Karyotyping
With minor adaptations, we followed the protocol of Sahara et al. (1999) for chromosome preparations. Specimens were anesthetized with carbon dioxide prior to dissection, and three legs were snap-frozen in liquid nitrogen and stored at À70°C for use in flow cytometry. Specimens were dissected in physiological solution containing 154 mmol/L NaCl, 5 mmol/L KCl, 1.8 mmol/L CaCl 2 , and 2.9 mmol/L NaHCO 3 . Gonads were removed and placed into hypotonic solution (0.075 mol/L KCl) for swelling (15-20 min), followed by fixation in Carnoy's solution (ethanol: chloroform: acetic acid = 6:3:1) for 20 to 30 min. Immediately after the removal of gonads, specimens were preserved in 99% ethanol and stored at À20°C for further use. Glass slides were incubated in acid ethanol (1% HCl: 96% ethanol = 1:100) for at least half an hour before use. Gonads were dispersed on drycleaned slides in a drop of acetic acid (60%) with tungsten needles and spread using a heating plate (40°C).
Chromosome preparations were first checked for metaphase spreads under a Nikon Eclipse E600 phase contrast microscope and then stained with 0.02 lg/mL DAPI in antifade solution based on DABCO (1,4-diazabicyclo[2.2.2] octane; Sigma-Aldrich). Pictures of chromosome spreads were taken with a Leica DFC 495 digital camera (Leica Microsystems, Wetzlar, Germany) mounted on a Leica DM 5000B fluorescence microscope. Single chromosomes were arranged according to their size using Adobe Creative Suite 4 software (Adobe Systems, San Jose, CA). Chromosome number was counted in at least three metaphase spreads per individual. Chromosome base number (x; i.e., the haploid number, which is occurring twice in a diploid and three times in a triploid) was calculated by dividing diploid chromosome number by two and triploid chromosome numbers by 3. Fundamental number of arms (FN; i.e., the total number of chromosome arms in a 2n metaphase complement) was determined for all species from representative karyotypes of as many sampling locations as possible.
In M. ticinensis, chromosome lengths (in lm) were measured in each one karyotype of nine asexual and five sexual individuals using the software MicroMeasure version 3.3 (Reeves 2000). We calculated mean chromosome length (MChL, equal to the sum of chromosome lengths divided by the number of chromosomes within one karyotype) and tested for significant difference in MChL among sexual and asexual individuals. Statistical tests were calculated in SigmaPlot v.12.5 (Systat Software Inc., San Jose, CA), applying a significance level of a = 0.05. When normality and equal variances were given, we applied parametric tests. Whenever one of these criteria failed, we used nonparametric tests instead (Reeves 2000).

Flow cytometry measurements
We applied flow cytometry to obtain relative genome-size estimates from leg muscle tissue using Bellis perennis (Asteraceae; 2C = 3.38 pg, Sch€ onswetter et al. 2007) as internal standard (i.e., measured together with each sample). Samples were treated following Suda et al. (2007). In detail, single legs and approx. 0.5 cm 2 of a B. perennis leaf were chopped with 500 lL of ice cold Otto 1 Buffer (0.1 mol/L citric acid, 0.5% Tween20 (Merck KGaA, Darmstadt, Germany)) in a small Petri dish. The suspension was then filtered through a 42-lm nylon mesh and incubated for a few minutes. Finally, 1 mL of staining solution (DAPI (4 lL/mL) and 2-mercaptoethanol (2 lL/ mL) in Otto 2 buffer (0.4 mol/L Na 2 HPO 4 •12 H 2 O)) was added. Fluorescence intensity of 3000 cell nuclei per sample was measured using a Partec CyFlow â space flow cytometer (Partec GmbH, M€ unster, Germany). Gating and peak analysis were performed automatically using the Partec Flo Max software. In a few samples that were defrosted too fast, fluorescence peaks showed increased noise on the left side. In these cases, left gates were adjusted manually to prevent overestimation of peak width.
Genome size measurements based on DAPI fluorescence have been criticized because of DAPI's AT-specific DNA-binding properties. Potential differences in genomic AT/GC content among species may introduce some bias to fluorescence intensity and thus genome size estimates. However, Johnston et al. (1999) showed that genome size estimates based on DAPI largely agree with estimates based on other fluorochromes. As we used a plant (B. perennis) as internal standard, we report relative genome size (i.e., fluorescence intensity of sample divided by fluorescence intensity of internal standard), instead of calculating nuclear DNA content in picograms. Relative genome size values were divided by ploidy level to obtain monoploid genome size (Cx), which is comparable among different ploidy levels. Fig. 4 was drawn in Sigmaplot v.12.5. Fig. 5 was generated using the R package ggplot2 (http://cran.r-project.org/web/packages/ggplot2/ index.html).

Results
For each species, an overview of sampling locations is given in Table 1, including number of female and male individuals sampled, mean chromosome number and relative genome size, ploidy level, elevation, and FN.

Mitochondrial phylogeny
Our phylogeny included eight target species, as well as four additional species mentioned in the specimen collection  section. All species formed monophyletic clusters supported by posterior probabilities above 0.95 (Fig. 2). In contrast, some deeper nodes were only weakly supported. The three sexual species formed a monophyletic group, while the parthenogenetic species did not. Mean net p-distances between species ranged between 3.00% and 16.22%, while p-distances within-species were equal or lower than 1.00% (Table 2). Mean within-group distances were significantly higher in sexual than in parthenogenetic species (Student's t-test: t = À3.369; df = 6; P = 0.015). Intraspecific splits occurred in M. helleri, M. hrabei, M. ticinensis, and M. tirolensis (Fig. 2). In M. ticinensis and M. tirolensis, these splits were congruent with reproductive mode and ploidy level, respectively. In M. helleri, the two resulting clades were partly congruent with differing chromosome number (Fig. 2 and Fig. 3). In M. hrabei, the two clades corresponded to sampling locations from Austria (KNB, VIE) and the Czech Republic (BRN; see Table 1).

Chromosome numbers
We found predominantly mitotic chromosomes in female gonads and meiotic chromosomes in male gonads, except in two male specimens (each one M. helleri and M. lehnhoferi), where both meiotic and mitotic chromosomes were found. In most sampling locations of the three sexual species (M. helleri, M. hrabei, and M. lehnhoferi), we consistently counted 2n = 52 chromosomes in females and n = 26 bivalents in males. The assembled karyotypes of these species revealed that they differed in the numbers of metacentric, submetacentric, and acrocentric chromosomes, and hence also in the fundamental number of arms (FN; Table 1, Fig. S1). In three locations of M. helleri, however, we found only 2n = 50 in females and n = 25 bivalents in males (HWD and RAX), and 2n = 51 chromosomes and n = 25 bivalents in one male from GIE ( Fig. 3; Table 1). Unfortunately, no unambiguously countable chromosome spread could be produced in female M. helleri specimens from GIE and RAX. All M. helleri specimens with 52 chromosomes fell into one well-supported mitochondrial subclade, including some specimens from locations RAX (all) and HWD (six individuals), which had only 50 chromosomes. The remaining individuals from HWD (five individuals) and GIE corresponded to the basal subclade within M. helleri (Fig. 3), having 50 chromosomes. In asexual species, chromosome number ranged from 2n = 46 to 2n = 3x = 78 (Fig. 2, Table 1). In M. fuscistylis, we detected intraspecific variation with 2n = 54 chromosomes across four locations (OBG, FOT, MAR, and SRK) but consistently 2n = 56 chromosomes in specimens from HIN ( Fig. 2 and Fig. S1). In all M. pallida specimens, we counted 2n = 3x = 78 chromosomes. In M. tirolensis, most sampling locations harbored exclusively diploid specimens (2n = 50 chromosomes), but at two locations (LDL and LDC), only triploids with 2n = 3x = 75 chromosomes were found (Fig. S1). We consistently found 2n = 50 chromosomes in specimens of M. engiadina.
In M. ticinensis, both males and females were found at location MIR, while exclusively females were found at the other five sampling locations (see Table 1), possibly indicating a pattern of geographic parthenogenesis. In the sexual population, specimens had 2n = 46 chromosomes plus a varying number of supernumerary elements, that is, B chromosomes, representing 2n = 46 (N = 3), 2n = 46 + 1B (N = 1), 2n = 46 + 4B (N = 3), and 2n = 46 + 5B (N = 1) in females, and n = 23 (N = 2) bivalents in male meiotic spreads (online supplementary material Fig. S1 and Tab. S1). In contrast, parthenogenetic specimens consistently had 2n = 46 chromosomes and no B chromosomes. In specimens from the sexual population, mean chromosome size was significantly larger compared with parthenogenetic specimens (Mann-Whitney U-test; U = 0, P = 0.003, N sex = 5, x sex = 4.26 AE 0.38 lm, N parth = 9, x parth = 3.33 AE 0.12 lm). Table 2. Genetic P-distances (in %) within (AESD) and among species. Rows 1-3 are sexual species, rows 4-8 (light gray background) are parthenogenetic species, and rows 9-12 (dark gray background) are species that were included in the phylogeny but not in other analyses due to low sample size.  Chromosome base numbers (x) significantly increased with elevation across all species (Pearson's r: r = 0.503, P < 0.001; Fig. 4A). However, sexual species did not contribute to this trend as no correlation was found among populations of sexual species alone (Pearson's r: r = À0.024, P = 0.919; Fig. 4A). In contrast, when only populations from parthenogenetic species were considered, the correlation was strongest (Pearson's r: r = 0.785, P < 0.001; Fig. 4A). No such correlation was found within any of the single species (data not shown).

Discussion
In this study, we investigated the link between parthenogenesis and polyploidy in Alpine representatives of the jumping-bristletail genus Machilis by generating relative genome-size estimates and karyotypes from geographically representative samples of three sexual and five presumably parthenogenetic species. By mapping relative genome size and chromosome numbers onto a mitochondrial phylogeny, we tackled the questions raised in the introduction and identified promising trajectories for future research.

Multiple origin of parthenogenesis
As this is the first phylogenetic framework in the genus Machilis, no 'a priori' hypotheses concerning species relationships and evolution of parthenogenesis were available.
Assuming that a reversal from parthenogenetic to sexual reproduction is unlikely following Dollo's law (Gould 1970), the nonmonophyly of parthenogenetic species indicated that asexuality originated at least five times independently (see Fig. 2). Multiple evolution of parthenogenesis within a genus has been demonstrated in other animal groups such as stick insects (Schwander and Crespi 2009) and lizards (Manriquez-Moran et al. 2014), and these groups offer exciting opportunities for studying evolutionary scenarios leading to parthenogenesis. In Alpine Machilis species, multiple origins of parthenogenesis may have been facilitated by Pleistocene glaciation cycles, which are thought to have increased the incidence of asexual reproduction in plants and animals (H€ orandl 2009). In fact, persistence on so-called nunataks (i.e., rocky peaks towering above the ice shield) has been suggested in two alpine Machilis species (M. fuscistylis : Janetschek 1956;M. pallida: Wachter et al. 2012) and may be linked to the onset of parthenogenesis in these species. In contrast, inner-Alpine persistence of glaciations at lower elevations is unlikely. In species occurring at lower elevations (e.g., M. engiadina, M. ticinensis, and M. tirolensis), parthenogenesis is thus more likely to have evolved during recolonization of the inner-Alpine area following persistence in peripheral refugia. However, in the absence of a robust time calibration, dating speciation events and origins of parthenogenesis remain elusive. In line with the concept of geographic parthenogenesis (Vandel 1928;H€ orandl 2009), parthenogenetic species included in our study mainly occurred in central parts of the Alps, whereas sexual species were mostly scattered close to the margin of the Alps and in surrounding lowland areas.

Incidence of polyploidy
Polyploidy was found in two parthenogenetic species. In M. tirolensis, populations with both, 50 and 75 chromosomes were identified. Because relative genome size in the latter was on average 1.58 times higher than in the former, the latter's triploid status is strongly supported. The second case of triploidy is hypothesized in M. pallida (2n = 3x = 78) based on the 1.5 ratio in chromosome number compared with sexual species. However, the direct sexual ancestor of M. pallida is unknown and possibly extinct, and thus, accurate inference of a phylogenetically unbiased genome size ratio is not possible. Still, a ratio of 1.34 was calculated using the average relative genome size across sexual species included in this study. Deviation from the expected 1.5 ratio might be explained by phylogenetic distance between M. pallida and these sexual taxa.
Since the other three parthenogenetic species were diploid, polyploidy is not always coupled with parthenogenesis in Machilis. Similarly, parthenogenetic diploids have been found in numerous other arthropod taxa, for example, brine shrimps (Maccari et al. 2013), weevils (Stenberg et al. 2003), and stick insects (Schwander and Crespi 2009). Our results are thus in line with the hypothesis, that in animals, polyploidy emerges secondarily, after the establishment of diploid parthenogenetic populations (Suomaleinen et al. 1987).

Opportunities for studying genome evolution in the genus Machilis
In this study, we found substantial variation in chromosome number and relative genome size among and within Machilis species. In the following, we highlight three peculiarities in our data that emphasize the potential of the genus Machilis as a study system for investigating the role of chromosomal rearrangements and genome size alterations in evolutionary diversification.

Genome downsizing
Monoploid genome size (Cx) was negatively correlated with chromosome base number (x) among species included in this study, and this correlation was strongest among parthenogenetic species. Moreover, while Cx decreased with elevation, chromosome base number tended to increase with elevation, and this effect was also stronger among parthenogenetic species. We cautiously interpret this as a signal of genome downsizing along an altitudinal gradient that is possibly stronger in parthenogenetic than in sexual species. Downsizing of nuclear DNA content is a known phenomenon in polyploid plants, where higher ploidy levels tend to have smaller C values than diploids of a given species (Leitch and Bennett 2004). Here, in contrast, we hypothesize that genome size and chromosome number might be shaped by environmental conditions and that parthenogenetic Machilis species might be more prone to this process than sexual congeners.
To date, specific mechanisms involved in genome size reductions are largely unclear. However, Hessen et al. (2010) proposed that reductions in genome size in eukaryotes may be explained by the need for streamlined genomes in environments selecting for rapid growth due to low food availability. As in alpine environments (i.e., above 2000 m a.s.l.) activity periods are short and food availability is low, we suggest that harsh climatic conditions selected for smaller genomes in alpine Machilis species. Alternatively, our interpretation might be misled by the fact that no exclusively alpine sexual species are included in our study. One reason for this is the unclear taxonomic situation in the genus Machilis , with several species being poorly delimited and thus not available for evolutionary research. Also, most putative Machilis species found in alpine habitats are presumably parthenogenetic, and one could expect that asexual reproduction itself is a crucial property enabling persistence in this habitat. Our results uncover a potential alternative hypothesis by pointing at smaller genomes (and more and thus smaller chromosomes) in alpine species. Future work should aim at validating this correlation in sexual alpine Machilis species to disentangle the relative contribution of genome size reduction versus asexual reproduction to the persistence in alpine environments.

Intraspecific chromosomal rearrangements
We found intraspecific variation in chromosome number within two species. In the parthenogenetic M. fuscistylis, one population (HIN) harbored two additional chromosomes compared with the other populations. As the fundamental number of arms differed between these subgroups, we excluded the possibility of chromosome fusion or fission, as both types of rearrangement would not alter FN. This assumption is further supported by relative genome size measurements, which confirmed a significantly higher 2C value in the HIN population. We propose that the duplication of one chromosome pair in HIN (Fig. S1) is more parsimonious than the parallel loss of one chromosome pair in the other populations. Also, as both karyotypic groups shared the same CO1 haplotype, the novel karyotype must have originated recently. We stress the importance of chromosomal rearrangements and karyotypic variation as a potential source of divergence among parthenogenetic lineages, which is in line with postulates by others (Sunnucks et al. 1998;Blackman et al. 2000).
Differing chromosome numbers in the sexual species M. helleri (52 vs. 50 chromosomes, one male with 51) showed that chromosomal rearrangements are not restricted to parthenogenetic Machilis species. As the fundamental number of arms was congruent across different karyotypes, the change in chromosome number is best explained by fission of one or fusion of two chromosomes, while the single male with 51 chromosomes may represent a hybrid karyotype or a chromosomal aberration. Alternatively, an additional chromosome in males compared with females could, for example, also indicate an XX/XYY' sex determination system. The geographical spacing of populations with 50 and 52 chromosomes (Fig. 3) is contrasted by an increase of 10% in relative genome size along a west-east gradient, resulting in larger genomes in populations having fewer chromosomes. Assuming that the ancestral chromosome number was 52, an increase in relative genome size might be explained as a by-product of the fusion of two chromosomes. Alternatively, this gradient may represent an independent phenomenon, possibly linked to an ecological cline or varying retrotransposon activity.
We hypothesize that the change in chromosome number in M. helleri could have built up a reproductive barrier, thus initiating genetic divergence. The evolutionary implications differ from those in parthenogenetic species (e.g., M. fuscistylis), where asexuality itself is a barrier to gene flow. This becomes evident in our mitochondrial phylogeny, where karyotypically differing individuals did not form reciprocally monophyletic groups but showed some incongruence in the populations RAX and HWD (Fig. 3), possibly a consequence of limited, ongoing gene flow between the karyotypic groups. Recently, it has been shown that chromosomal rearrangements are key factors promoting speciation during secondary contact in several species pairs of rodents (Castiglia 2014). Similarly, the pattern seen in M. helleri could be interpreted as an early stage of speciation triggered by a chromosomal rearrangement.

Occurrence of B chromosomes
Varying numbers of B chromosomes were found in the sexual population of M. ticinensis. Interestingly, relative genome size significantly increased with the number of B chromosomes within this population, but the significant difference in relative genome size between sexual and parthenogenetic populations could not be explained by their presence alone. This is consistent with the results of Trivers et al. (2004), who found approx. 35% larger genomes in plants with reported B chromosomes compared with plants without B chromosomes. Because B chromosomes are thought to include a high amount of transposable elements (Camacho 2005), repeated transposition of mobile elements to the original chromosome set and subsequent proliferation could explain increased nuclear DNA content in the sexual population of M. ticinensis. Studying B chromosomes is an excellent opportunity to learn more about the dynamics of evolving genomes (Houben et al. 2014), and thus, M. ticinensis represents a prime study system for future research targeting the impact of B chromosomes on genome-size evolution.