# Intraspecific Phenotypic Variation and Morphological Divergence of Strains of Folsomia candida (Willem) (Collembola: Isotomidae), the "Standard" Test Springtaill

Thomas Tully, Mikhail Potapov

► To cite this version:

Thomas Tully, Mikhail Potapov. Intraspecific Phenotypic Variation and Morphological Divergence of Strains of Folsomia candida (Willem) (Collembola: Isotomidae), the "Standard" Test Springtaill. PLoS ONE, 2015, 10 (9), pp.e0136047. 10.1371/journal.pone.0136047. hal-01271716

HAL Id: hal-01271716

https://hal.sorbonne-universite.fr/hal-01271716

Submitted on 9 Feb 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

RESEARCH ARTICLE

# Intraspecific Phenotypic Variation and Morphological Divergence of Strains of *Folsomia candida* (Willem) (Collembola: Isotomidae), the "Standard" Test Springtail

Thomas Tully<sup>1,2*</sup>, Mikhail Potapov<sup>3</sup>

1 Institute of ecology and environmental sciences—Paris (iEES Paris, UMR 7618), CNRS, UPMC Univ Paris 06, Sorbonne Universités, Paris, France, 2 ESPE de Paris, Université Paris 4, Sorbonne Universités, Paris, France, 3 Department of Zoology and Ecology, Moscow State Pedagogical University, Moscow, Russia

* thomas.tully@upmc.fr

## Abstract

We describe and compare the external morphology of eleven clonal strains and one sexual lineage of the globally distributed *Folsomia candida*, known as "standard" test Collembola. Of the 18 morphological characters studied, we measured 14 to have significant between-strains genetic variations, 9 of these had high heritabilities (>78%). The quantified morphological polymorphism was used to analyse the within-species relationships between strains by using both a parsimony analysis and a distance tree. These two detailed morphological phylogenies have revealed that the parthenogenetic strains grouped themselves into two major clades. However the exact position of the sexual strain remains unclear and further analysis is needed to confirm its exact relationship with the parthenogenetic ones. The two morphologically based clades were found to be the same as the ones previously described using molecular analysis. This shows that despite large within-strain variations, morphological characters can be used to differentiate some strains that have diverged within a single morphospecies. We discuss the potential evolutionary interpretations and consequences of these different levels of phenotypic variability.

## Introduction

Springtails belong to a very ancient group (Class Collembola, ~410 Mya) [1] which has colonised nearly all terrestrial habitats [2]. This taxon comprises about 8000 described species and is considered to be moderately diversified [2,3]. New species are described regularly, mostly on the basis of external morphology [4]. But the real number of existing species of Collembola on Earth is difficult to predict although a ballpark estimation of about 50,000 species has been proposed [5]. It seems also reasonable to assume that the current species' diversity is considerably underestimated, especially if one takes into account molecular studies that have recently revealed the existence of extensive cryptic genetic diversity in several widely distributed species or complexes of species [3,6–11].

<figure><img src="image_3.png" /></figure>

OPEN ACCESS

**Citation:** Tully T, Potapov M (2015) Intraspecific Phenotypic Variation and Morphological Divergence of Strains of *Folsomia candida* (Willem) (Collembola: Isotomidae), the "Standard" Test Springtail. PLoS ONE 10(9): e0136047. doi:10.1371/journal.pone.0136047

**Editor:** Peter Shaw, Roehampton university, UNITED KINGDOM

**Received:** April 24, 2015

**Accepted:** July 29, 2015

**Published:** September 10, 2015

**Copyright:** © 2015 Tully, Potapov. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** Data are available from Dryad (http://dx.doi.org/10.5061/dryad.64m57).

**Funding:** This study received financial support of the Ministry of Education and Science of the Russian Federation (Project No 6.632.2014/K) and of the Angence Nationale de la Recherche ANR-09-PEXT-011 and from UPMC and CNRS (France).

**Competing Interests:** The authors have declared that no competing interests exist.
To recognise, describe and define species in this group, taxonomists rely essentially on the chaetotaxy [4]. Using the chaetotaxy, they describe and compare the number, position and shape of external sensillar elements, searching for diagnostic characters even in minute species [11].

Intraspecific character variation is commonly observed in Collembola chaetotaxic studies [12–24] and for pragmatic reason it is normally considered as problematic rather than informative [25]. Some common setae vary considerably in number between individuals or even from one side to the other on the same individual. They also “float” making the determination more complicated [26]. Other setae are much more stable and are useful for interspecific identification [4,25,27].

The parthenogenesis is widely known in Collembola [28,29] and the morphological species concept which is used for practical purposes in the springtails is especially relevant for such a mode of reproduction for which the biological species concept can not be applied. In parthenogenetic species, the gene flow is disrupted between lineages. Each strain thus forms an independent evolutionary unit, as do different species that have recently diverged [30,31]. Different strains may have diverged owing to genetic drift of different selective pressures if they evolve in contrasted environments. Thus studying the infraspecific morphological character variation in a parthenogenetic species may be of interest for studying the early steps of morphological divergence during an on-going speciation event and understand species and their origins in asexual organisms [31].

Our work focuses on the evolution of infraspecific character variation in the parthenogenetic springtail *Folsomia candida*, one of the most commonly used laboratory Collembola [32]. This springtail is cosmopolitan [9,32] and its range spans every terrestrial ecozones. In Europe, it occupies mostly protected soils (flower pots, greenhouses, caves) though it can be collected in many open habitats being more frequent in more disturbed sites. This morphospecies is known to reproduce mainly asexually [32] although some populations are sexual [33].

Due to its worldwide distribution and to the ease with which it can be raised and maintained in the laboratory, this species became a standard micro-arthropod and is now used in many laboratories as a model organism in different domains of soil biology [32]. It is widely used especially as a model organism for ecotoxicological testing purposes [34–37].

Several strains from different geographical origins have been collected and used by various laboratories [9]. But most experiments are conducted on one [35,38–40] or a few clonal strains [41–43]. If the genetic homogeneity of different individuals belonging to a single strain is a clear benefit for standardisation in the laboratory experiments, it can also lead to a delusion as regards the relevance of comparisons that can be made between studies made on different strains. In other words, because of disclosed cryptic genetic diversity in *F. candida*, this Collembola may not be an ideal “standard” as it is supposed to be [44].

Indeed, some infraspecific genetic polymorphism has been described using allozyme techniques, RAPD markers or gene sequences [9,41,42]. In a previous study where 11 strains were compared, we found that this within-species genetic diversity is organised in two major lineages [9]. Although some genetic differences are found between different strains within each lineage, major and consistent differences have been observed between these two clades. It was later shown that two contrasted biodemographic strategies evolved along with the early divergence of these two major branches of the evolutionary tree. These two main strategies differ in their overall reproductive and survival potentials and on their level of adaptive plasticity [45,46]. One can easily see that making use of strains belonging to one branch or the other of this species evolutionary tree may strongly influence the outcome of any laboratory experiment.

Our aim here is to study within this commonly used morphospecies (*F. candida* sensu stricto, see below) whether external morphology can be used to make sense of the cryptic genetic diversity and less cryptic life history diversity, which have previously been described.
More precisely we aim to address the following four points:

i. **The description of morphological variation.** How stable is the external morphology in this morphospecies? Can we see some morphological trait variation and if so what are the traits that are found to vary? Does this variation put in question the taxonomical status of this morphospecies?

ii. **The quantification of morphological variation.** By comparing the morphology of individuals belonging to different strains, we can describe the patterns of morphological variation within and amongst lineages of a single species. This will enable us to know whether the traits that vary between strains also vary between individuals within a lineage. In other words, it will enable us to partition the total morphological variances of traits into their additive genetic components (between-lineage variance) and their random environmental or developmental components. This is required to quantify the heritabilities of the morphological traits that are found to vary.

iii. **The evolutionary significance of morphological variation.** Estimating the heritabilities of traits is important since it permits one to estimate the level of morphological divergence between the strains and to know how far these strains are on the early stages of an ongoing speciation process. Moreover, the quantification and partitioning of the morphological variation can be used to build up a phylogeny of the different lineages, and to question whether this phylogeny is congruent with the ones derived from molecular analysis. This will show us if the morphology, which is supposed to carry some reliable information for taxonomy, is also relevant for disentangling the evolutionary relationships of different lineages within a single morphospecies.

iv. **The practical use of morphological variation.** The present work also aims at searching for characters that could be used to reliably identify some lineages of this commonly used species. This could enable laboratories to easily verify whether they are working on common clades or on different ones.

# Material and Methods

## Species' taxonomical description

The springtail *Folsomia candida* was first described in 1902 by Willem [47] from a single individual found in a Belgian cave. It was later morphologically characterised in greater detail [48]. Afterwards this species was repeatedly redescribed following the understanding of Stach (for details see [49]). So far it is not considered to be a good morphospecies-Potapov et Gao [50] proposed to use two terms *F. candida* sensu lato and *F. candida* sensu stricto. The former is more tentative as including all forms or species closely related to "typical" *F. candida* greatly differing in body size, number of chaetae on furca, foil-chaetae and other features. *F. candida* sensu lato is more known from Asia and America. *F. candida* s.str. was defined basing on European strains kept in the laboratory of the first author of the paper, its diagnosis is much narrower [50]. In the present study we deal with *F. candida* s.str.

## Lineages examined

Twelve laboratory strains of *F. candida* have been examined, and labelled AP, BR, BV, DK, GB, GM, HA, PB, SH, TO, US, WI. All of them except SH reproduce asexually and have been kept in the laboratory since 1999. Apart from the strain SH, all were previously studied according to
their geographical origin and molecular phylogeny [9] and to some of their life-history traits [45,46]. As mentioned above, these previous studies have shown that these strains can be grouped into two major clades. The first clade (A) encompasses five lineages (AP, GB, HA, BR & BV) with a low reproductive potential but a relatively high longevity. The second clade (B) is composed of six lineages (DK, GM, PB, TO, US & WI) with a higher reproductive potential and a higher mean death rate [45,46]. As regards the strain SH, it comes from a population in Shanghai (China) and is composed of male and female individuals.

To exclude potential age-dependent morphological variability, only young adults were used for this study. They were all grown from eggs, received the same food and maintained in similar containers at the same temperature (21°C) to standardise conditions, which could influence the morphology. Note that SH was not cultivated in the same laboratory and we cannot exclude the fact that slightly different conditions could influence their morphology. All of the studied SH were females. For each strain, about twelve young adults having approximately the same size (~1.2/1.8 mm) were taken from these standard cultures and stored in alcohol.

A total of 66 individuals of *F. candida* with about 8 to 12 individuals per strain have been examined. More abundant material could give more characters to be used to separate the strains statistically.

Note that for each strain, the collected individuals were separated into two tubes and the 2*12 = 24 tubes were anonymised by the first author and then given to the second author for examination. The strain identity was later reattributed to the 24 groups of individuals after their character description was completed.

## Morphological characters

The springtails were mounted on microscope slides in Phoera liquid (200 g of chloral hydrate, 20 g of glycerol and 30 g of gum arabic dissolved and mixed into 50 g of distilled water) and their morphology was carefully studied. The specimens were first examined to verify that their chaetotaxy corresponded to the diagnostic characters of the morphospecies [49]. Then 18 external morphological characters easy to observe with a standard light microscope on individuals mounted on a slide were systematically examined. These characters are listed and described in Table 1 and the main discriminant characters are illustrated in Fig 1. Almost all of these characters have been somewhat used in species' taxonomy of the genus *Folsomia* [27], but their taxonomical values are different and depend on group of species. Additional sensilla on dorsal side of Ant.III (not shown) was not used because of its high variability within a strain. In addition, it is often difficult to recognise and calculate seta-like sensilla on this side of Ant. III. Many other characters were left out of our consideration here because they predictably displayed no or very low variability between our strains since we dealt with one species. We tried to observe every character on each studied individual although all the characters were not always visible. Obvious aberrations were not taken into consideration.

## Statistical analysis

**Character analysis.** We analysed the variability of each character using either a binomial linear model for dichotomous characters such as the presence of an additional post-labial setae (character n°2 in Table 1) or a Poisson model for count characters with large number of values such as the number of setae on Ant.I (n°3) which varies for instance between 18 and 32. But most count characters displayed an almost binomial distribution. By that we mean that for these characters, most of the observed variation spreads between two values with some rare extra-values. For instance, about 95% of the observed springtails had either 3 or 4 setae on upper subcoxeae of leg 2 (n°6) whereas only 4% had two and one individual had five setae. For
<table><caption>**Table 1. Characters description and level of variation.** The whole range of character variation is given. Data in parenthesis represent variants occurring singly. The whole dataset is available in S1 Dataset.</caption><thead><tr><th>Body part</th><th>(Character number), name and description.</th><th>Range of character variation (outliers); Threshold (/)</th><th>Type of model and overall differences between clones and</th><th>Estimation of the character heritability h<sup>2</sup>.</th></tr></thead><tbody><tr><td rowspan="5">Head</td><td>(01) The second guard setae along posterior edge of Post-Antennal Organ (PAO). The second seta is in group of guard setae of PAO. It is located near the middle part of PAO and is normally longer than the others.</td><td>At a level with other guard setae / behind other guard setae</td><td>Binomial (1 for "at the level"), χ<sup>2</sup><sub>11</sub> = 133, P &lt; .001</td><td>94%</td></tr><tr><td>(02) Nb. of post-labial setae on one side. All strains bear at least 4 large basal post-labial setae on one side. The variant of 5 setae means that an additional smaller seta appears between the first and second basal setae.</td><td>4–5; 4 / 5</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 135, P &lt; .001</td><td>78%</td></tr><tr><td>(03) Nb. of setae on Ant.I. Counted on one side. Only common setae are counted, sensilla and basal microsensilla are out of consideration. Basal microsensilla are shorter and thinner than common setae and are located dorsally and ventrally on the basal part of the segment.</td><td>18–32; ≤ 24 / &gt;24</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 27, P = .005</td><td>4.4%</td></tr><tr><td>(04) Nb. of sensilla on Ant.I. Total numbers of sensilla are reported here. There is invariably one long sensilla associated with generally 3 (2 to 4) short sensilla. Most of the character variation being between 4 and 5, we considered to have only two states: 5 or less than 5.</td><td>4–5; 4 (rarely 5) / 5 in all individuals of strain.</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 79, P &lt; .001</td><td>86%</td></tr><tr><td>(05) Sensilla on lateral side of Ant.III. Nb. of sensilla (0, 1, or rarely 2) located just posterior to lateral outer sensillum of antennal organ of third antennal segment. On right and left sides (usually visible on one side only). The additional sensillum is located just posterior to lateral outer sensillum of antennal organ of third antennal segment. It is shorter and thinner than the latter one.</td><td>0–1 (2); absent in all individuals of strain / present in all or some individuals of strain</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 107, P &lt; .001</td><td>94%</td></tr><tr><td rowspan="6">Thorax</td><td>(06) Nb. of setae on upper subcoxa of Leg II</td><td>3–4</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 19.4, P = .054</td><td></td></tr><tr><td>(07) Nb. of setae on lower subcoxa of Leg II</td><td>6–11</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 5.8, P = .88</td><td></td></tr><tr><td>(08) Nb. of setae on upper subcoxa of Leg III</td><td>4–9</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 5.6, P = .90</td><td></td></tr><tr><td>(09) Nb. of setae on lower subcoxa of Leg III</td><td>7–14</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 17, P = .08</td><td></td></tr><tr><td>(10) Tooth on claw. Inner tooth is visible in laterally positioned claw and is normally estimated on the third pair of leg. The state "present" means that the tooth is well visible; irregularity in curvature of inner edge of claw is considered the absence.</td><td>absent or indistinct / present</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 137, P &lt; .001</td><td>88%</td></tr><tr><td>(11) Nb. of ventral setae on Th.III. Counted on one side.</td><td>1–4; 1 on both or one side / 2 or more on both or one side</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 113, P &lt; .001</td><td>78%</td></tr><tr><td rowspan="5">Abdomen</td><td>(12) Nb. of anterior setae on dens.</td><td>22–38</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 24.2, P &lt; .001</td><td>6.7%</td></tr><tr><td>(13) Nb. of posterior setae in basal group on dens. All setae on posterior side of dens can be divided to a basal and two proximal groups: outer and inner. The basal group is located on strictly posterior part of dens lacking crenulations. As a rule, the mentioned groups are separated by gaps; in some individuals they are less distinct.</td><td>4–9; ≤ 7 / ≥ 6</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 37, P &lt; .001</td><td>12%</td></tr><tr><td>(14) Nb. of posterior setae in outer-proximal group on dens.</td><td>2–4</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 88, P &lt; .001</td><td>80%</td></tr><tr><td>(15) Nb. of posterior setae in inner-proximal group on dens.</td><td>2–4</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 20, P = .04</td><td>15%</td></tr><tr><td>(16) Nb. of apical setae on posterior side of manubrium.</td><td>(0–1) 2–4 (5); ≤ 2 / ≥ 2</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 57, P &lt; .001</td><td>94%</td></tr><tr><td>(17) Nb. of posterior (dorsal) setae in lateral parts of manubrium (right/left).</td><td>(0) 1–2 (3)</td><td>Binomial, χ<sup>2</sup><sub>11</sub> = 122, P &lt; .001</td><td>85%</td></tr><tr><td>(18) Nb. of anterior setae on manubrium. Measured on both sides. One or two setae on lateral sides of manubrium are not included to the group of anterior setae.</td><td>11–27; ≤ 18 / ≥ 18</td><td>Poisson, χ<sup>2</sup><sub>11</sub> = 59, P &lt; .001</td><td>12%</td></tr></tbody></table>

doi:10.1371/journal.pone.0136047.t001

this kind of character, we grouped together the outliers and rare values with their nearest most common ones in order to analyse the data as a binomial variable. In the former case, 2 and 3 are grouped together and 4 and 5 form the second group. This was done for characters n°4, 5, 6, 10, 11, 14, 15, 16 and 17. But in order to show the range of variation of the observed characters we have conserved the different observed possible values in the graphs (Figs 2, 3 and 4).
<figure><img src="image_4.png" /><figcaption>**Fig 1. External morphology.** External morphology of *Folsomia candida* displaying the main observed characters and their variation. The character numbers correspond to the one in Table 1. *bms*, basal microsensilla.</figcaption></figure>

doi:10.1371/journal.pone.0136047.g001
Head

<figure><img src="image_5.png" /><figcaption>**Fig 2. Head.** Head morphological character state distribution (Table 1) for the studied strains of *F. candida*. The strains are aligned as follows: the new SH sexual strain is plotted on the left hand side. Then, based on previous work, we have grouped together the clones AP to HA since they are known to be genetically affiliated. The right hand side of the graphs represents in alphabetic order the clones DK to WI that are known to belong to another distinct clade. For each strain the raw observations are plotted with little vibration so that the points do not overlap. For lateral characters, the left side characters are figured as circles and the right side ones as squares. The estimated means and 95% confidence intervals (see main text for details, especially for characters treated as binomial variables) are plotted for each strain. The horizontal grey dotted line underlines the overall mean trait value.</figcaption></figure>

doi:10.1371/journal.pone.0136047.g002

The outliers in these graphs are plotted close to their nearest more dominant values so that the predicted mean and confidence interval derived from the binomial model could be plotted on top of the raw data (Figs 2, 3 and 4).

Binomial or Poisson linear models were used to search and test for any potential differences between strains. We used the *glm* (generalised linear model) function from the R statistical package to fit these models [51]. For each character, an appropriate *glm* model was fitted to the data with strain labels as covariate to test whether the trait differs globally between the strains
<figure><img src="image_6.png" /><figcaption>**Fig 3. Thorax.** Thorax morphological character state distribution for the studied strains of *F. candida*. See Fig 2 caption and Table 1 for details.</figcaption></figure>

doi:10.1371/journal.pone.0136047.g003

(p values in Table 1). This analysis was used to sort out the characters that displayed significant differences between strains from those with no overall differences.

Similar *glm* models with strain labels as covariate but without the intercept were used to derive the predicted means and 95% confidence intervals for each character and strain (Figs 2, 3 and 4, Table 1). Only the characters that displayed some significant genetic variation were later used to measure their heritability and to build up a distance cladogram and a phylogenetic tree.

**Characters' heritabilities.** We estimated the broad sense heritabilities of the characters that were found to differ significantly between strains. For clonal organisms, the relevant measure of genetic variance is the broad-sense heritability defined as the ratio of the among-clone component of variance to the total phenotypic variance: $H^2 = \sigma_G^2/\sigma_T^2$ where $\sigma_G^2$ is the genetic
Abdomen

<figure><img src="image_7.png" /><figcaption>**Fig 4. Abdomen.** Abdomen morphological character state distribution for the studied strains of *F. candida*. See Fig 2 caption and Table 1 for details.</figcaption></figure>

doi:10.1371/journal.pone.0136047.g004
A) Strict-consensus cladogram

<figure><img src="image_8.png" /></figure>

Mean values of the scaled and
centered parameter values

B) Cluster dendrogram

(UPGMA-like distance tree)

<figure><img src="image_9.png" /></figure>

C) Character numbers

**Fig 5. Phylogeny.** A) Strict consensus tree obtained from a parsimony analysis of tree based on the binomilisation of the scaled and centred morphological data (black and white colours of the trait values in C). B) Cluster dendrogram made on the scaled measurements of the morphological traits that have been shown to vary between strains using average method. The approximately unbiased (au) p-values (expressed as %) computed by multiscale bootstrap resampling indicate for each node how strong the cluster is supported by data.

doi:10.1371/journal.pone.0136047.g005

variance (between strain variance) and σ²_T is the total trait variance [52]. We used the *glmmPQL* function from the *MASS* package to do such estimations. *glmmPQL* fits a generalised linear mixed model using penalised quasi-likelihood. In our case this model has been fitted to binomial or Poisson data depending on the characters.

We used no fixed effect but introduced in the model the strain label variable as a random factor to partition the overall character variance between “between-strain” variance and “within-strain” variance which enables the model to estimate the proportion of total character variation due to differences between strains, thus providing a surrogate of the broad sense heritability of the characters [52,53] (Table 1).

**Phylogenetic analysis.** We first studied the relationship between the strains using a cladistic approach. The character coding in the matrix of character was done as follows. We scaled and centred the matrix of the estimated mean per strain values of the morphological variables and transformed this matrix into a binomial matrix using the sign of each transformed measurement. For each character, the strains whose means are lower than the average mean have a state “0” whereas the other ones have a state “1” (black and white colours in Fig 5C). Here, we only used the morphological characters that were found to significantly vary between the strains. This matrix was transformed into an object of class *phyDat* (*phangorn* library in R) and then processed by the *bab* function (*phangorn* library) which finds all the most parsimonious trees after examining all possible trees [54,55]. We then used the *consensus* (*ape* library) function in order to obtain the strict consensus tree from the most parsimonious trees. This consensus tree was then artificially rooted between the two previously described clades (clades labelled A and B in Fig 5A) for graphical homogeneity [9].

**Cluster analysis.** The cladistic cladogram was complemented with a phenetic distance tree dendrogram. To build-up this dendrogram from our morphological data, we used a
hierarchical cluster analysis (*hclust* function in program R 2.15) with an average linkage method which is the one used for UPGMA (Unweighted Pair Group Method of Aggregation). This was performed on the matrix of the scaled (centred and standardised) mean values for each strain of the characters that were shown to vary significantly (Fig 5C). Multiscale boot-strap resampling (*pvclust* function from the *pvclust* library) was used to provide approximately unbiased p-values (expressed as %) for each node, which indicate how strong the clusters are supported by the data. This enabled to highlight the clusters that are highly supported by the data [56] (Fig 5B).

# Results

## Characters' variability and heritability

All individuals could be identified as the previously recognized *F. candida* morphospecies sensu strictu [50].

Some phenotypic variation was observed for all 18 characters under study but the level of genetic variation varied a lot between characters.

No genetic variation was found for four characters: the number of setae on all subcoxa of second and third legs (n° 6–9, Fig 3) strongly varied depending on the specimen but this variation was not found to be structured between the strains.

Some significant genetic variation associated with low heritability ($h^2<15\%$) was found for five characters (n° 3, 12, 13, 15, 18, Table 1), four of which were counts of setae (Poisson distributed data). These low heritabilities seemed to be due to the large within-strain character variation with large overlapping of the character distribution between the clones. For the number of anterior setae on dens (n°12), all the clonal strains were found to be very similar, the only observed genetic variation coming from the sexual strain SH which had on average less setae than the other strains. Aberrant asymmetry of dens often occurs–many individuals bear a different number of setae on the right and left dens (24 and 30, 27 and 35, as examples). Some differences between strains could probably be found if more specimens were studied.

Nine characters have high heritabilities (n° 1, 2, 4, 5, 10, 11, 14, 16, 17, $h^2>78\%$, Table 1), all of them being binomial characters. But even for these characters, some within-strain variation could be observed with only few exceptions such as for instance the alignment of a guard setae on PAO for SH, BR, GB, TO and US (N°1), the presence of additional sensilla on Ant. 3 (n° 5) for SH, AP, BR, BV, GB, WI, the presence of a clawtooth (n° 10) for SH, BR, DK or the number of ventral setae on thorax 3 for strains SH and DK (n° 11). Because of this within-strain variation, the assignment of a unique and clear-cut value for each clone-character pair is not straightforward.

## Phylogenetic and distance trees

The parsimony analysis yielded two primary trees. We knew from a previous study that the eleven clonal strains group themselves into two distinct clades [9]. After artificially rooting our morphologically based phylogenetic tree, we found that our strains organised themselves again into two clades that correspond to the previously described clades: AP, GB, HA, BR and BV in clade "A" and DK, GM, PB, TO, US & WI for clade "B" (Fig 5A). The strict consensus tree has two multifurcations on clade B. In this analysis the sexual strain SH branches itself close to the clone GB. Based on morphology, the clone AP branches at the base of clade A. This was also the case in our previous molecular phylogeny [9].

The overall topology of our cluster dendrogram (Fig 5B) was found to be very similar to the phylogenetic cladogram: the two lineages (A & B) are clearly supported. And even within these two lineages, some similarity can be observed: the similarity between clones PB and GM or
between BR and GB is supported for instance by both approaches (Fig 5). The main difference between the two topologies is the position of the sexual strain SH which roots itself outside the two clades in the dendrogram whereas it is closely linked to clone GB in the cladogram. This discrepancy is probably due to the character n°12 which provides no information in the cladistic analysis whereas it moves SH away from the other strains on the distance tree (Fig 5C).

## Discussion

### Taxonomical status of strains

All the strains under study fit quite well the modern taxonomical understanding of the species *F. candida* sensu strictu [50]. All strains have characteristic long furca, the same sensillar chaetotaxy, ventral setae on third thoracic segment, 13-29 setae on anterior side of manubrium and other characters that define well this species (for more details see the redescriptions mentioned above and [50]). In our individuals, the following characters differ a little from that is written in the diagnosis compiled by Potapov for this species in more traditional and wider understanding [27]: number of apical setae on posterior part of manubrium (1-5 versus 2), number of anterior setae on manubrium (11-27 versus 16-32), sensilla on first antennal segment (4-5 versus 3), ventral setae on third thoracic segment (1-4+1-4 versus 2-3+2-3). No doubt, our strains don't represent all possible variability of this species sensu lato.

### Variability and heritabilities of the characters

All characters vary strongly depending on individuals. This yields to considerable overlapping between clones in almost all characters. But most of the characters studied here displayed some level of genetic variability. The range of within-strain variation of the four characters with no genetic variation was quite high and this may explain the absence of any significant heritability. The heritabilities of the other characters varied from comparatively low values to very high values (from ~7% up to 94%). It is interesting to note that the heritabilities never reached 100%, which means that all of the heritable morphological characters that we have studied were not absolutely stable for every lineage. To our knowledge these are the first measurements of morphological characters' heritabilities in a Collembola.

The large range of measured heritabilities can reflect the fact that the evolution and divergence of morphological characters among isolated clonal strains happens at different rates: some characters have already diverged a lot while for others the emergence of genetic variability may still be at its beginning if it is ever going to evolve. But it can also result from different levels of phenotypic canalisation among the characters. For the same level of genetic variability, a high canalisation is associated with a higher heritability due to a lower level of within-lineage character variation. Overall, this character diversification reflects the polychaetotic or oligo-chaetotic evolutionary tendencies within this group of *F. candida* strains (apart from characters n°1 and 10).

### Strain identification by morphology

According to our analysis, between 9 and 13 of 18 morphological characters can be used to discriminate the strains depending on the pair of strains one wants to compare. These characters may be strongly age-dependent and should only be used among the individuals of the same age or instar [57]. Most of them are a number of setae on different parts of the body.

These morphological characters could be profitably used to verify that two or more strains belong to the same clade or to quantify their morphological difference. This could be useful, for instance, for assessing the validity of comparison of multiple ecotoxicological tests performed
on different strains. It could complement the use of molecular tools for genetic characterisation of the strains [41]. These morphological characters could also be used if one wants to run laboratory experiments where several strains are mixed together in competition experiments for example [58]. Using an appropriate character and pair of clones, it is possible to identify a strain using a single individual. For instance, one individual of the strains BR, BV and GB can be clearly distinguished from an individual of strain US by counting fewer or more than 24 setae on Ant.I (n° 3) respectively.

But for most comparisons, drawing a robust conclusion means studying more than one individual. To distinguish two strains the choice of characters will depend on how many individuals can be measured. If one wants to identify individuals belonging to several strains mixed together in a single population, one will either need to use very discriminant characters or to breed the isolated individuals to make more measurements on less discriminant characters on their numerous descendants.

## Patterns of morphological and molecular evolution

Using analysis based on the genetically variable morphological characters, we found that the variability of morphological characters carries some phylogenetic signal and can be used to study the patterns of evolution of the different lineages within a single species. The two methods we used to build-up the trees (cladistic and distance) produced very similar results.

The studied strains first split into two distinct clades and then diverged within each clade. We have again found the two major clades that have been previously found using molecular phylogeny [9].

The congruence between the topology of the distance tree and the cladistic tree is striking. The topology is similar at a broad level (the two main clades) but also at a closer look. In our previous work using some genetic markers [9], the genetic signal was insufficient to provide us with some details regarding the fine topology of the tree within each clade. More work needs to be done using molecular analysis to tell whether the congruence between the tree based on molecular and morphological data persists at a thinner phylogenetic scale. But we can conclude from this work that morphological characters, even if they display a lot of variation within a population or within a strain, can provide in this species reliable and meaningful information regarding the identification of clonal strains and the reconstruction of their evolutionary relationships.

## The evolution of setosity

If we compare on average the evolutionary tendency of the chaetotaxy of the strains between the two clades, it appears that clade A tends to be oligochaetotic while the second clade tends to be more polychaetotic. This pattern is striking for characters n° 2, 3, 5, 16, 18 while the reverse tendency only appears for character n°11. The lineages with a slow life history (A) are less hairy than the fast ones (B). Is this morphological evolution adaptive and if yes what is the adaptive significance of this pattern? The two clades have very contrasted life-history strategies and are thus believed to have adapted to contrasted environments [45,59].

During this evolutionary divergence, it is possible that multiple small-effect mutations controlling the level of setosity on different parts on the body have arisen independently and have been selected differently within the two major clades. Under this adaptive scenario, the level of setosity could be considered as an adaptation, and selection should also act to reduce the developmental variance by favouring robust genotypes producing more canalised phenotypes whose level of setosity is adapted to a specific environment [60]. But from this viewpoint, the adaptive significance of this morphological evolution remains to be determined. However the observed
pattern of phenotypic variation could also result from the fixation of a single or few mutations of relatively large effects that control the level of setosity globally, or at least on several parts of the body as has been shown to be the case for *Drosophila* [61,62].

In such a case, the observed difference in the mean level of setosity between the two clades would not necessarily be adaptive since it could simply result from the fixation by chance of one of few alleles controlling the number of setae on different parts of the body in the common ancestor of one of the two clades.

Lastly one can question the genetic basis on which the evolution of setosity in *Folsomia candida* depends. Much work has been done on the genetic mechanisms underlying bristle patterns in Diptera and especially in *Drosophila* [63,64]. The *achaete-scute* (ac-sc) gene codes for a transcription factor which plays a central role in determining the development of the fly's macrochaetes [63]. The position of bristles on the thorax is determined by the spatial expression of this *achaete-scute* gene and this expression is itself controlled by many independent cis-regulatroy elements [64–66]. Mutations in these cis-regulatory regions can affect some of the thoracic bristles in *Drosophila* [61,66,67] and the variation of the number and arrangement of these bristles in different species of Diptera seems to be driven by the evolution of these cis-regulatory sequences [63,64].

One can make use of these results when studying the genetic basis of the evolution of setos-
ity in *Folsomia candida*. Genomic tools developed for springtails and the recently sequenced
and analysed transcriptome of *F. candida* [68] could be used in order to search for an *ac-sc*
homologue and its associated cis-regulatory regions which must exist in *Folsomia* genome
given that both the structure and function of the *ac-sc* genes are known to be conserved in
many animal species [63,64,69]. This could then be used for studying whether the observed
phenotypic divergence in the setosity rely, as for *Drosophila*, on mutations in cis-regulatory ele-
ments regulating the transcription of a collembolan homologous of the *ac-sc* gene [64].

More work is required (i) to study the selective advantages that could be provided by various
levels of setosity under different environmental conditions, (ii) to study the adaptive signifi-
cance of the observed morphological development noise and (iii) to search for the genetic
determinants underlying these phenotypic variations and assess whether these morphological
changes are driven by the evolution of cis-regulation.

A speciation in progress?

We knew that *F. candida* morphospecies comprises a great number of life-history diverging evolutionary lineages [45,59,70]. We observed some phenotypic variation on many morphological traits and found that the lineages also diverged morphologically. This is a good evidence that morphological speciation is in progress.

We observed that the evolution of new phenotypes in diverging lineages of *F. candida* is
accompanied by a high developmental noise (non genetic phenotypic variability) since for each
lineage, the diverging morphological traits are still unstable and are partly overlapping. How
this variability may evolve during what could be an on-going speciation remains to be
addressed. Moreover, this divergent process may be blurred and slowed down by genetic trans-
fers between the diverging lineages driven by the sexual strains if they can intercross with some
of these lineages. But the precise position and role of these sexual strains in the evolutionary
history of this species remains at present unsettled.

Cryptic diversity and potentialities of morphology in Collembola

A high level of cryptic diversity within already known Collembola is often postulated in the
current publications [7,8,71,72]. At whole for Collembola, the unexpected genetic divergence
was revealed within homogeneous morphospecies. As a result, several common and widely distributed species passed to the “complexes” of many potential species: *Parisotoma notabilis*, *Friesea grisea*, *Heteromurus major*, *Ceratophysella denticulata* and others. It was concluded that modern morphological taxonomy underestimates the number of species. We showed in our study that the morphological characters could be successfully found and used to identify even the intraspecific forms, which are the strains in our study. The considerable overlapping between clones in almost all characters makes the identification of clones less exact but possible. With this aim in mind we must render the morphological noise as small as possible, particularly by culturing the individuals in the same conditions and by examining individuals of the same age, that can partly be applied to the species mentioned above.

## Supporting Information

### S1 Dataset.

(ZIP)

## Author Contributions

Conceived and designed the experiments: TT MP. Performed the experiments: TT MP. Analyzed the data: TT. Contributed reagents/materials/analysis tools: TT MP. Wrote the paper: TT MP.

## References

1. Whalley P, Jarzembowski EA. A new assessment of *Rhyniella*, the earliest known insect, from the Devonian of Rhynie, Scotland. Nature. 1981, May 28; 291: 317.

2. Hopkin SP. Biology of the Springtails (Insecta: Collembola). Oxford: Oxford University Press; 1997.

3. Cicconardi F, Fanciulli PP, Emerson BC. Collembola, the biological species concept and the underestimation of global species richness. Mol Ecol. 2013, Oct 1; 22 (21): 5382–96. doi: 10.1111/mec.12472 PMID: 24112308

4. Deharveng L. Recent advances in Collembola systematics. Pedobiologia. 2004; 48: 415–33.

5. Hopkin SP. Collembola: the most abundant insects on earth. Antenna. 1998; 22: 117–21.

6. Cicconardi F, Nardi F, Emerson BC, Frati F, Fanciulli PP. Deep phylogeographic divisions and long-term persistence of forest invertebrates (Hexapoda: Collembola) in the North-Western Mediterranean basin. Mol Ecol. 2010, Jan; 19 (2): 386–400. doi: 10.1111/j.1365-294X.2009.04457.x PMID: 20015142

7. Porco D, Potapov M, Bedos A, Busmachiu G, Weiner WM, Hamra-Kroua S, Deharveng L. Cryptic diversity in the ubiquist species *Parisotoma notabilis* (Collembola, Isotomidae): A long-used chimeric species? PLoS One. 2012, Sep 26; 7 (9): e46056. doi: 10.1371/journal.pone.0046056 PMID: 23049931

8. Zhang F, Yu D, Luo Y, Ho SYW, Wang B, Zhu C. Cryptic diversity, diversification and vicariance in two species complexes of *Tomocerus* (Collembola, Tomoceridae) from China. Zoologica Scripta. 2014, May 13; 43 (4): 393–404.

9. Tully T, D'Haese CA, Richard M, Ferriere R. Two major evolutionary lineages revealed by molecular phylogeny in the parthenogenetic Collembola species *Folsomia candida*. Pedobiologia. 2006; 50 (2): 95–104.

10. Schneider C, Cruaud C, D'Haeese CA. Unexpected diversity in Neelipleona revealed by molecular phylogeny approach (Hexapoda, Collembola). Soil Organisms. 2011; 83: 383–98.

11. Schneider C, D'Haese CA. Morphological and molecular insights on *Megalothorax*: the largest Neelipleona genus revisited (Collembola). Invertebrate Systematics. 2013; 27 (3): 317–64.

12. Yoshii R. On *Lobella sauteri* and its variability. Ann. spel. Inst. Japan. 1994; 12: 1–9.

13. de Mendonça C, Arlé R. Variação quetotáxica em *Brachystomella agrosa* Wray, 1953 (Collembola, Neanuridae). Boletim do Museu nacional Nova Série Zoologia. 1992 (351): 1–9.

14. Wetton MN. Morphological variation in British *Folsomia quadrioculata* Tullberg (Collembola: Isotomidae): a multivariate study. Systematic entomology. 1987; 12 (2): 257–70.

15. Mari Mutt JA. New records for *Corynothrix borealis* and a study of its variability (Collembola: Entomobryidae: Orchesellinae). International Journal of Entomology. 1984; 26 (4): 369–77.
16. Pitkin BR. Variation in some British material of the *Onychiurus armatus* group (Collembola). Systematic Entomology. 1980; 5 (4): 405-26.

17. Lawrence PN. Variation in the Palaeotropical Collembola species *Cryptopygus fasciatus* with three new synonyms and notes on the subgenus *Rhodanella*. *Revue d'écologie et de biologie du sol*. 1978; 15 (2): 265-71.

18. Goto HE, Ögel S. Variation in the mucro of *Folsomia* (Collembola: Isotomidae). The Entomologist. 1961; 94: 105-7.

19. Lawrence PN. Variation in *Hypogastrura* (H.) *viatica* (Tullberg, 1872) (Collembola). The entomologist. 1961, Jun: 151-4.

20. Rapoport EH, Cutini de Izarra DE. On the variability of *Hypogastrura manubrialis* (Tullb.)(Collembola). Annals and MAgazine of Natural History. 1962; 5 (52): 205-8.

21. Rusek J. Zup Variabilität einiger morphologischen Merkmale bei Onychiurus pseudostachianus gisin, 1956 (Collembola, apterygota). Věstník Československé společnosti zoologické (Acta Societatis Zoologicae Bohemicae). 1964; 28 (1): 44-7.

22. Christiansen K, Culver D. Geographical variation and evolution in *Pseudosinella hirsuta*. Evolution. 1968; 22 (2): 237-55.

23. Gough HJ. On variation in *Arrhopalites caecus* (Tullberg)(Collembola: Sminthuridae). Entomologist's Monthly Magazine. 1973; 108: 205-9.

24. Prabhoo NR. Morphological variations in the local populations of the soil collembola *Proisotoma (Claviosotoma) fitchioides* (Denis 1947)(Insecta) in Kerala. Entomon. 1976; 1 (1): 67-70.

25. Grow A, Christiansen K. Chaetotaxy in *Folsomia* (Collembola Isotomidae) with special reference to nearctic species. *Revue d'Ecologie et de Biologie du Sol*. 1976; 13 (4): 611-27.

26. Grow AB, Christiansen K. Chaetotaxy in nearctic *Friesea* (Collembola Neanurinae) with notes on taxonomic use of chaetotaxy. *Revue d'Écologie et de Biologie du Sol*. 1974; 11 (3): 377-96.

27. Potapov MB. The use of macrochaetotaxy in taxonomy of palaearctic *Folsomia* (Collembola, Isotomidae). Pedobiologia. 2000; 44 (3-4): 234-239.

28. Petersen H. Sex ratios and the extent of parthenogenetic reproduction in some collembolan populations. In: Dallai R, editors. *First International Seminar on Apterygota*. Accademia delle Scienze di Siena detta de'Fisiocritici, Siena; 1978. p. 19-35.

29. Chernova M, Potapov B, Savenkova YY, Bokova I. Ecological significance of parthenogenesis in Collembola. Entomological Review. 2010, Mar; 90 (1): 23-38.

30. Lost sex: the evolutionary biology of parthenogenesis. Schön I, Martens K, van Dijk P, editors. Dordrecht, Heidelberg, London, New York: Springer; 2009.

31. Birky CW Jr, Barraclough TG. Asexual speciation. In: Schön I, Martens K, van Dijk P, editors. Lost sex: the evolutionary biology of parthenogenesis. Dordrecht, Heidelberg, London, New York: Springer; 2009. p. 201-16.

32. Fountain MT, Hopkin SP. *Folsomia candida* (Collembola): a "standard" soil arthropod. Annual Review of Entomology. 2005, Jan; 50: 201-22. PMID: 15355236

33. Krogh PH. Toxicity testing with the collembolans *Folsomia fimetaria* and *Folsomia candida* and the results of a ringtest. Miljøstyrelsen; 2009.

34. Bur B, Crouau C, Bianco B, Gandois G, Probst P. Toxicity of Pb and of Pb/Cd combination on the springtail *Folsomia candida* in natural soils: Reproduction, growth and bioaccumulation as indicators. Science of The Total Environment. 2012; 414 (0): 187-197.

35. Crouau Y, Chenon P, Gisclard C. The use of *Folsomia candida* (Collembola, Isotomidae) for the bioassay of xenobiotic substances and soil pollutants. Applied Soil Ecology. 1999, May; 12 (2): 103-11.

36. Environnement Canada. Méthode d'essai biologique: essai de mesure de la survie et de la croissance de collemboles exposés à des contaminants dans le sol. Ottawa, Canada: Centre des sciences et technologies environnementales. Direction générale des sciences et de la technologie Environnement Canada; 2007.

37. ISO 11267. Soil quality—Inhibition of reproduction of Collembola (*Folsomia candida*) by soil contaminants. Genève: Inter Stand Org; 2014.

38. Waagner D, Heckmann LH, Malmendal A, Nielsen NC, Holmstrup M, Bayley M. Hsp70 expression and metabolite composition in response to short-term thermal changes in *Folsomia candida* (Collembola). Comparative Biochemistry and Physiology A—Molecular & Integrative Physiology. 2010, Oct; 157 (2): 177-83.

39. de Boer TE, Birlutiu A, Bochdanovits Z, Timmermans MJTN, Dijkstra TMH, Van Straalen NM, et al. Transcriptional plasticity of a soil arthropod across different ecological conditions. Molecular Ecology. 2011, Mar; 20 (6): 1144-54. doi: 10.1111/j.1365-294X.2010.04985.x PMID: 21255171
40. Pfeffer SP, Khalili H, Filser J. Food choice and reproductive success of *Folsomia candida* feeding on copper-contaminated mycelium of the soil fungus Alternaria alternata. *Pedobiologia*. 2010; 54 (1): 19-23.

41. Diogo JB, Natal-da-Luz T, Sousa JP, Vogt C, Nowak C. Tolerance of genetically characterized *Folso-mia candida* strains to phenmedipham exposure—A comparison between reproduction and avoidance tests. *Journal of Soils and Sediments*. 2007; 7 (6): 388-92.

42. Chenon P, Rousset A, Crouau Y. Genetic polymorphism in nine clones of a parthenogenetic collembolan used in ecotoxicological testing. *Applied Soil Ecology*. 2000, Apr; 14 (2): 103-10.

43. Stam EM, de Leemkule MAV, Ernsting G. Trade-offs in the life history and energy budget of the parthenogenetic Collembolan *Folsomia candida* (Willem). *Oecologia*. 1996; 107 (3): 283-92.

44. Crouau Y, Cazes L. What causes variability in the *Folsomia candida* reproduction test? *Applied Soil Ecology*. 2003; 22 (2): 175-80.

45. Tully T, Ferrière R. Reproductive flexibility: genetic variation, genetic costs and long-term evolution in a collembola. *PLoS One*. 2008; 3 (9): e3207. doi: 10.1371/journal.pone.0003207 PMID: 18791644

46. Tully T, Lambert A. The evolution of postreproductive life span as an insurance against indeterminacy. *Evolution*. 2011, Oct; 65 (10): 3013-20. doi: 10.1111/j.1558-5646.2011.01347.x PMID: 21967441

47. Willem V. Note préliminaire sur les Collemboles des Grottes de Han et de Rochefort. *Ann. Soc. Ent. Belg*. 1902; 46: 275-83.

48. Stach J. The Apterygotan fauna of Poland in relation to the world fauna of this group of Insects family: Isotomidae. Nakl. Polskiej Akademii Umiejętności; 1947.

49. Potapov MB. Synopses on Palaearctic Collembola, Volume 3: Isotomidae. Staatliches Museum für Naturkunde; 2001.

50. Potapov M, Yan G. *Folsomia* of China I—fimetaria group (Collembola: Isotomidae). *Annales de la Société Entomologique de France*. 2012; 48 (1-2): 51-6.

51. Ihaka R, Gentleman R. R: a language for data analysis and graphics. *Journal of computational and graphical statistics*. 1996; 5 (3): 299-314.

52. Lynch Michael, Walsh Bruce. Genetics and analysis of quantitative traits. Sunderland, Mass.: Sinauer; 1998.

53. Bolker BM, Brooks ME, Clark CJ, Geange SW, Poulsen JR, Stevens MH, White JS. Generalized linear mixed models: a practical guide for ecology and evolution. *Trends Ecol Evol*. 2009, Mar; 24 (3): 127-35. doi: 10.1016/j.tree.2008.10.008 PMID: 19185386

54. Paradis E. Analysis of phylogenetics and evolution with R. New York: Springer; 2012.

55. Schliep KP. phangorn: phylogenetic analysis in R. *Bioinformatics*. 2011, Feb 15; 27 (4): 592-3. doi: 10.1093/bioinformatics/btq706 PMID: 21169378

56. Suzuki R, Shimodaira H. Pvclust: an R package for assessing the uncertainty in hierarchical clustering. *Bioinformatics*. 2006, Jun 15; 22 (12): 1540-2. PMID: 16595560

57. Potapov MB. Age-dependent variability of ventral chaetotaxy of manubrium in *Folsomia spinosa* Ksene-man, 1936. Fauna i ecologia nogokhvostok, Moscow, Nauka. 1984: 26-8.

58. Le Bourlot V, Tully T, Claessen D. Interference versus exploitative competition in the regulation of size-structured populations. Am Nat. 2014; 184 (5): 609-23. doi: 10.1086/678083 PMID: 25325745

59. Tully T, Lambert A. The evolution of postreproductive life span as an insurance against indeterminacy. *Evolution*. 2011, Oct; 65 (10): 3013-20. doi: 10.1111/j.1558-5646.2011.01347.x PMID: 21967441

60. Flatt T. The evolutionary genetics of canalization. The Quarterly Review of Biology. 2005, Sep; 80 (3): 287-316. PMID: 16250465

61. Stern DL, Orgogozo V. Is genetic evolution predictable? *Science*. 2009; 323 (5915): 746. doi: 10.1126/science.1158997 PMID: 19197055

62. McGregor AP, Orgogozo V, Delon I, Zanet J, Srinivasan DG, Payre F, Stern DL. Morphological evolution through multiple cis-regulatory mutations at a single. *Nature*. 2007; 7153 (587-U6): doi: 10.1038/nature05988

63. Simpson P, Woehl R, Usui K. The development and evolution of bristle patterns in Diptera. Development. 1999; 126 (7): 1349-64. PMID: 10068629

64. Simpson P. The stars and stripes of animal bodies: evolution of regulatory elements mediating pigment and bristle patterns in *Drosophila*. *Trends Genet*. 2007, Jul; 23 (7): 350-8. PMID: 17499383

65. Gómez-Skarmeta JL, Campuzano S, Modolell J. Half a century of neural prepatterning: the story of a few bristles and many genes. *Nat Rev Neurosci*. 2003, Jul; 4 (7): 587-98. PMID: 12838333

66. Marcellini S, Simpson P. Two or four bristles: functional evolution of an enhancer of scute in Drosophila. *PLoS Biol*. 2006; 4 (12): e386. PMID: 17105353
67. Gibert JM, Marcellini S, David JR, Schlötterer C, Simpson P. A major bristle QTL from a selected population of *Drosophila* uncovers the zinc-finger transcription factor poils-au-dos, a repressor of achaete-scute. Dev Biol. 2005, Dec 1; 288 (1): 194–205. PMID: 16216235

68. Faddeeva A, Studer RA, Kraaijeveld K, Sie D, Ylstra B, Mariën J, et al. Collembolan transcriptomes highlight molecular evolution of Hexapods and provide clues on the adaptation to terrestrial life. PLoS One. 2015; 10 (6): e0130600. doi: 10.1371/journal.pone.0130600 PMID: 26075903

69. Grens A, Mason E, Marsh JL, Bode HR. Evolutionary conservation of a cell fate specification gene: the Hydra achaete-scute homolog has proneural activity in *Drosophila*. Development. 1995, Dec; 121 (12): 4027–35. PMID: 8575303

70. Stam EM, de Leemkule MAV, Ernsting G. Trade-offs in the life history and energy budget of the parthenogenetic Collembolan *Folsomia candida* (Willem). Oecologia. 1996; 107 (3): 283–92.

71. Rougerie R, Decaëns T, Deharveng L, Porco D, James SW, Chang C-H, et al. DNA barcodes for soil animal taxonomy. Pesquisa Agropecuaria Brasileira. 2009; 44 (8): 789–802.

72. Torricelli G, Carapelli A, Convey P, Nardi F, Boore JL, Frati F. High divergence across the whole mitochondrial genome in the "pan-Antarctic" springtail *Friesea grisea*: evidence for cryptic species? Gene. 2010, Jan 1; 449 (1–2): 30–40. doi: 10.1016/j.gene.2009.09.006 PMID: 19782734