Mountain Refugia Play a Role in Soil Arthropod Speciation on Madagascar: A Case Study of the Endemic Giant Fire-Millipede Genus Aphistogoniulus

To elucidate the speciation mechanisms prevalent within hotspots of biodiversity, and the evolutionary processes behind the rise of their species-rich and endemic biota, we investigated the phylogeny of the giant fire-millipede genus Aphistogoniulus Silvestri, 1897, a Malagasy endemic. This study is the first comprehensive (molecular and morphological) phylogenetic study focusing on millipede (class Diplopoda) speciation on Madagascar. The morphological analysis is based on 35 morphological characters and incorporates ten described as well as two newly described species (A. rubrodorsalis n. sp. and A. jeekeli n. sp.) of Aphistogoniulus. The molecular analysis is based on both mitochondrial (COI and 16S), and nuclear genes (complete 18S rDNA), together comprised of 3031 base pairs, which were successfully sequenced for 31 individual specimens and eight species of Aphistogoniulus. In addition to the null-model (speciation by distance), two diversification models, mountain refugia and ecotone shift, were discovered to play a role in the speciation of soil arthropods on Madagascar. Mountain refugia were important in the speciation of the A. cowani clade, with three species occurring in the Andringitra and Ranomafana Mountains in the southeast (A. cowani), the Ambohijanahary and Ambohitantely Mountains in the mid-west (A. sanguineus), and the Marojejy Mountain in the northeast (A. rubrodorsalis n. sp.). An ecotone shift from the eastern rainforest to the unique subarid spiny forest of Mahavelo was discovered in the A. vampyrus - A. aridus species-pair. In the monophyletic A. diabolicus clade, evidence for divergent evolution of sexual morphology was detected: species with greatly enlarged gonopods are sister-taxa to species with normal sized gonopods. Among the large-bodied Spirobolida genera of Madagascar, Colossobolus and Sanguinobolus were found to be close sister-genera to Aphistogoniulus. Forest destruction has caused forest corridors between populations to disappear, which might limit the possible resolution of biogeographic analyses on Madagascar.


Introduction
Madagascar is one of the world's centers of endemism and is considered a biodiversity hotspot [1]; with its highly diverse ecosystems, it mirrors a small continent rather than an island. Although now located off the coast of eastern Africa, until 88 million years ago Madagascar formed a landmass with India [2]. The Madagascar-India landmass separated from continental Africa in the early Mesozoic era, 158-160 million years ago [3], although land bridges towards Antarctica might have existed up to a more recent date [4]. The long isolation of Madagascar from other continents resulted in mixed-origins of the Malagasy fauna: neo-endemics whose ancestors colonized the island and quickly diversified [5][6][7]; and truly ancient endemics, whose ancestors were present before Madagascar became an island [8][9][10][11].
The isolation of Madagascar makes the island an ideal model region for the study of diversification in species, because unlike on most continents, introgression or admixture events with lineages that evolved outside the area can be almost ruled out [12]. Recent studies, almost exclusively on vertebrates, present different hypotheses of diversification mechanisms on Madagascar (see overview in [12]), including mountain refugia [13,14], retreatdispersal watersheds [15], and river barriers [16]. Unfortunately, some speciation patterns might already be artificially modified due to habitat destructions by humans, especially in the transitional zones [17].
Currently, the study of speciation events on Madagascar focuses primarily on vertebrates (e. g. [14,[18][19][20][21]). Although terrestrial invertebrates represent the largest percentage of Madagascar's diversity, little is known about their biogeography. The existing studies of Madagascar's invertebrates focus on taxa more likely to be dispersers, either through their ability to fly [6,22,23] or to wind disperse [24]. Non-flying soil invertebrates are, with some notable exceptions in mollusks [25,26], generally ignored. For Madagascar, the general absence of studies of taxa with few dispersal abilities is problematic because, unlike in mobile taxa, certain dispersal mechanisms can be ruled out, which provides a clearer view of the past [27]. Many soil invertebrates are good model taxa for biogeographic studies because they possess poor dispersal abilities, are often confined to now discontinuous habitats, and show slow growth and low fecundity [28]. For example, in the Australian wet tropics non-flying invertebrates have a much greater number of subregional endemics compared to vertebrates (50% vs. 9%, [29]). Furthermore, invertebrates have smaller habitat requirements than vertebrates [30], and are more likely to have a continuous history of occupation on a site than a mammal species [27].
The incorporation of less vagile lineages into the studies of speciation events on Madagascar may be useful because they are phylogenetically older, can survive in smaller habitats during climatically unfavorable times, and might represent relictual elements, living witnesses of a climatically different past.
This study utilizes the giant fire-millipede genus Aphistogoniulus Silvestri, 1897 as a test case for invertebrate speciation on Madagascar. As in other groups of organisms [31], species representing the largest, most widespread and most colorful Malagasy millipede genus were described first. For the Spirobolida from Madagascar, this was the large-bodied fire-millipede genus Aphistogoniulus [32], with its first species being described almost 130 years ago [33]. Species of Aphistogoniulus feature a striking pitchblack/blood-red aposematic color, sometimes accompanied by golden legs, and can reach a body length of up to 185 mm. A complete taxonomic revision of Aphistogoniulus was recently undertaken [34], providing a foundation for the phylogenetic study undertaken here. The wide distribution of the genus renders Aphistogoniulus a suitable model taxon to test different methods of species diversification on Madagascar (as proposed by [12]). Species of the genus can be found in montane and lowland rainforests all the way from the southeastern tip of the island, expanding north, leaving out only the far north of Madagascar [34]. Two species almost never live in sympatry, but replace each other in adjacent areas. Only a single species, A. aridus Wesener, 2009, lives in the dry spiny forest.
Although a complete revision of the genus Aphistogoniulus was undertaken recently, two additional undescribed species were discovered in more or less recent collections, including the pet trade. These new discoveries prompted the attempts undertaken in this study to construct a phylogeny of all species of the genus Aphistogoniulus, in combination with all other large-bodied Spirobolida genera from Madagascar. This study is one of only a handful of molecular systematic studies of millipedes [35] and the first dealing with millipedes in the southern hemisphere. It is investigated which diversification mechanisms, such as mountain refugia [13] or retreat-dispersion watersheds [15], most-likely played a role in the speciation of the endemic rainforest soil biome on Madagascar.

Nomenclatural Acts
The electronic version of this document does not represent a published work according to the International Code of Zoological Nomenclature (ICZN), and hence the nomenclatural acts contained in the electronic version are not available under that Code from the electronic edition. Therefore, a separate edition of this document was produced by a method that assures numerous identical and durable copies, and those copies were simultaneously obtainable (from the publication date noted on the first page of this article) for the purpose of providing a public and permanent scientific record, in accordance with Article 8.1 of the Code. The separate print-only edition is available on request from PLoS by sending a request to PLoS ONE, Public Library of Science, 1160 Battery Street, Suite 100, San Francisco, CA 94111, USA along with a check for $10 (to cover printing and postage) payable to ''Public Library of Science''. In addition, this published work and the nomenclatural acts it contains have been registered in ZooBank, the proposed online registration system for the ICZN. The ZooBank LSIDs (Life Science Identifiers) can be resolved and the associated information viewed through any standard web browser by appending the LSID to the prefix ''http://zoobank. org/''. The LSID for this publication is: urn:lsid:zoobank.org:pub:BE53B87F-9FDA-4B75-B461-A7110D373A99; urn:lsid:zoobank.org:pub:BE53B87F-9FDA-4B75-B461-A7110D373A99.

Morphological Phylogenetic Analysis
Taxon selection. Because the closest relative of Aphistogoniulus is unknown, several genera were selected as outgroups.  Wesener, 2009. For the ingroup, all ten described species as well as two newly discovered species of Aphistogoniulus, distributed above the whole island of Madagascar ( Fig. 1) were added to the character matrix. All material was either the type series, or a specimen which was directly compared to the types. For all specimens studied see Table 1.
Character selection. A character matrix was created using Mesquite 2.6 ( [38], see Table 2). A total of 35 characters were selected, focusing on the male copulatory devices (gonopods, see also Supporting Information S5). Two character sets of Aphistogoniulus, the endotergum ( = underside of the body rings) and the female vulva (see Supporting Information S3) were illustrated for the first time (see Supporting Information S4). Due to a lack of sufficient material, endotergum characters could not be added to the matrix (see Supporting Information S6). All characters were scored as 'unordered' with equal weights. Drawings were produced using a standard camera lucida mounted on an Olympus SZH10 dissecting microscope. For scanning electron microscopy, samples were cleaned and dehydrated in an ethanol series and air-dried overnight before being sputter coated (Denton Vacuum Desk IV) for 180 seconds. SEM micrographs were taken using a Zeiss (Leo) EVO SEM, based at the FMNH.
Phylogenetic analysis. A maximum parsimony Branch and Bound tree search was conducted in PAUP* 4.0b10 [39]. Accelerated transformation was used as character optimization criterion. Six of the 35 characters were parsimony uninformative. The search yielded a single shortest tree with a length of 35. Furthermore, a bootstrap analysis, incorporating 2000 pseudoreplicates, was undertaken in PAUP 4.0b10 under the tree-bisection-reconnection (TBR) branch-swapping algorithm.  [34].

Molecular Phylogenetic Analysis
Specimen sampling and DNA extraction. Freshly collected specimens were euthanized through freezing, which relaxes the specimen and allows for further study without breaking. Two legs were removed immediately and stored in 95% ethanol which was changed once after approximately four weeks. The specimen was then transferred to 75% ethanol. The other 22 specimens were previously stored in 75% ethanol and were collected between 1995 and 2006 (see Supporting Information S1 and Supporting Information S8). From these specimens, muscle from the male gonopods was removed and transferred to 95% ethanol. All tissue samples came from male specimens, with the exception of M. maximus B and A. infernalis E which came from mature females. Total genomic DNA was extracted using the DNAeasyH Blood & Tissue kit from Qiagen TM , following the manufacturer's extraction protocol. The leftover DNA extraction product as well as the complete second elution will be stored at the Field Museum's cryogenic storage facility. To allow a later identification and verification of the species after a taxonomic update, voucher specimens were deposited (see Table 3).
DNA amplification and sequencing. Three genes, the mitochondrial partial 16S rRNA gene, the cytochrome oxidase c subunit I (COI) partial gene, and the complete nuclear 18S rRNA gene, were amplified using the polymerase chain reaction (PCR) [41]. All three genes and gene fragments were successfully used in previous millipede studies and provided good resolution up to family level [11]. Amplification reactions were carried out on a MJ Researchß PTC 200 thermal cycler. Temperature profile and summaries of the amplification reactions can be found in previous studies [11]. Negative and positive controls were included in every PCR setup. The PCR product was purified with a commercial kit (QIAquick ß PCR Purification Kit, Qiagen GmbH). Cycle sequencing for sequencing was conducted with BigDye on a Bionad DyadH DNA Engine, cleaned with ethanol on an Eppendorf ß centrifuge 5810R and sequenced on a Hitachi/AB 3730 automatic DNA sequencer, using the same primer sets as for PCR or, in the case of the 18S rDNA, specific primer pairs (see Supporting Information S9). Sequencing reads were assembled and proofread with the program Seqman TM II (DNASTAR, Inc.) while the identity of all new sequences was confirmed with BLAST searches [42]. All used primer sequences are listed in Supporting Information S9, while all new sequences were deposited in GenBank (see Table 3 for accession numbers).
Maximum Parsimony phylogenetic analysis. The program PAUP*4.0b10 [39] was used for maximum parsimony (MP) analyses using the TBR branch swapping algorithm with an unenforced 'MaxTrees' option. The number of parsimony informative characters in the dataset was 645. Total number of shortest length trees with 2114 steps was 15, which all showed identical topologies on the species-level. A strict consensus tree was built out of all 15 shortest-length trees. To test if no tree island was overlooked during the heuristic search, a Winclada-Asado ß Ratchet analysis [44] was conducted under Asado version 1.7 ß . The Ratchet parameters were set to: number of iterations = 200 (default), number of trees to hold at each iteration = 5 (default = 1), number of characters to sample = 350 (default = 303). All other settings were within the default parameters of Asado 1.7.
To assess statistical support for hypothesized clades, 1,000 bootstrap pseudoreplicates were calculated with the program PAUP*4.0b10, under the TBR branch swapping algorithm, with one tree held at each step and without enforcing maximum trees.
Bayesian phylogenetic analysis. Bayesian analyses were conducted using MrBayes v3.1.2 [45]. Appropriate DNAsubstitution models were determined separately for all analysed three genes under the Bayesian information (BIC) [46] implemented in jModeltest 0.1.1 ( [47,48]. In the Bayesian analysis, the combined data set was partitioned (see MrBayes manual for details) to allow unlinked models and model parameters for all three genes. Following the MrBayes manual, the most appropriate substitution model was specified without fixing the parameter values, allowing them to vary during the annealing process. This method leads to more conservative results but more realistic posterior probabilities. The Bayesian analysis was conducted by computing 3,000,000 Monte Carlo Markov chain (MCMC) generations in two parallel runs each with three cold chains and one hot chain. Trees were sampled every 100 generations. The number of burn-in generations was determined by manual inspection of the likelihood score over the 3,000,000 generations.
Maximum Likelihood phylogenetic analysis. Using RAxML 7.0.4 [49] a fast maximum likelihood (ML) method including bootstrap analyses with 10,000 pseudoreplicates was performed. The dataset was partitioned to allow different parameters for each of the three analysed genes under GTRMIX parameter settings (see RAxML 7.0.4 manual for details).

Analysis of the Historical Biogeography of Aphistogoniulus
In order to explain the current biogeographic distribution of the analyzed eight species of Aphistogoniulus we used RASP (formerly S-DIVA) to reconstruct the habitat ranges of the ancestral taxa [50]. Habitat classifications were defined based on the collection site of the studied species (see Fig. 1), following the atlas of the vegetation of Madagascar. Habitat A refers to seasonal dry forests, habitat B refers to montane rain forests, habitat C refers to lowland rainforests and habitat D refers to mid elevation rain forests. It is important to note that the outgroup taxon Colossobolus semicyclus was chosen primarily for rooting the phylogeny, not for determining the ancestral habitat of the genus Aphistogoniulus. A Bayesian binary MCMC analysis was run using default settings: number of chains = 10, frequent of sample = 100, discard samples = 100, temperature = 0.1, maximum number of areas = 4, state frequencies = fixed (Jukes-Cantor), and an across-site rate variation = equal. However, we increased the number of cycles from 50,000 to 100,000. All results in detail are listed in Supporting Information S11.

Morphologic Phylogenetic Analysis
The first split in the tree occurs between the members of the tribe Pachybolini (Epibolus pulchripes, Madabolus maximus) and the rest of the large-bodied Spirobolida from Madagascar (Fig. 2). Characters supporting this split are found in the anatomy of the head (Supporting Information S5, c1) and mouthparts (c2), as well as the special shape of the vulva (c3, c4, c5), and is also supported by a bootstrap value (100%) in the analysis (Fig. 2).
The next well-supported split divides the southern Malagasy genus Corallobolus from the three genera Sanguinobolus, Colossobolus and Aphistogoniulus (Fig. 2: clade A). All members of clade A share a strongly elongated coxite (Supporting Information S5, c19) and a disc-shaped telopodite of the posterior gonopods (S5, c20). Sanguinobolus & Colossobolus form the sister-group to Aphistogoniulus (Fig. 2). Aphistogoniulus is monophyletic in the tree (Fig. 2, clade B) and supported by apomorphies such as the elongation of the coxite process of the anterior gonopods (S5, c15) and the location of the telopodite retrorse process on the anterior gonopods (S5, c16).
Three clades are present within Aphistogoniulus, the A. cowani, A. erythrocephalus and A. diabolicus clades. The genus is divided by a basal split between the A. cowani and A. erythrocephalus sister clades ( Fig. 2: clade C) and a clade comprising A. corallipes and the A. diabolicus clade (Fig. 2: clade F). Clade (C) is supported by the presence of a retrorse projection on the anterior side of the main branch of the posterior gonopods (S5, c25). The main apomorphy of the A. cowani clade is the lateral membranous fringe which is developed into a strong process on the main branch of the posterior gonopods (S5, c24). The A. erythrocephalus clade is corroborated by a sulcate fringe present on the apex of the posterior gonopod main branch (Supporting Information S2).
Three characters support the A diabolicus clade: the special shape (S5, c29) and the sharp, pointed ending (c33) of the basal branch of the posterior gonopod, as well as the presence of a apico-mesal fringe on main branch of the posterior gonopod (c30). This clade also includes the only species of Aphistogoniulus adapted to arid environments, A. aridus, which is the sister taxon of A. vampyrus (Fig. 2). The current position of A. corallipes is not well-supported in the tree, being in a basal position towards the A. diabolicus clade.

Molecular Phylogenetic Analyses
The molecular phylogenetic tree obtained by the parsimony analysis (Fig. 3) closely resembles the morphological one (Fig. 2). Aphistogoniulus and Sanguinobolus form a well-supported group (Fig. 3,  clade A), distinct from the other Malagasy genera Spiromimus (Spiromiminae) and Madabolus (Pachybolini). The monophyly of Aphistogoniulus (Fig. 3, clade B) is also well supported with a bootstrap value of 99%. Inside Aphistogoniulus, we observe a basal split between the species of the A. cowaniand A. erythrocephalus group (clade C) and the A. diabolicus clade (Fig. 3, clade F). The A. cowani clade (Fig. 3, clade D) with A. cowani and A. sanguineus is strongly supported (99%) and clearly distinct from the clade incorporating A. erythrocephalus and A. hova (Fig. 3, clade D). A. sanguineus and A. cowani differ moderately in their mitochondrial DNA (COI: 14-16%, 16S: 4.4-5.7%, see Supporting Information S10). The A. diabolicus clade (Fig. 3, clade F) is statistically wellsupported (83%) and also incorporates A. corallipes. All species of the clade (A. jeekeli, A. corallipes, A. infernalis and A. vampyrus) come out as monophyletic with strong statistical support. The northernmost species of the group, A jeekeli n. sp., is the sister-taxon to all others, followed by A. corallipes, which is the sister taxon to the A. vampyrus/A. infernalis species-pair. The genetic difference between A vampyrus and A. infernalis is moderate (COI: 15-17%, 16S: 5.7-6.1%), despite their quite distinct copulation legs, especially the posterior gonopods [34]. In A. infernalis the specimens from the littoral rainforest of Sainte Luce are the sister group (differences COI: 6.6%, 16S: 2.0%) to populations from the isolated rainforest at Grand Lavasoa and from the Andohahela Mountain chain (Fig. 3).
For the Bayesian analysis, jModeltest indicated the following models as the appropriate nucleotide substitution models for the three analysed genes, including the following parameters: 18S rDNA -Selected model: TPM3+G; 16S rDNA -Selected model: TrN+G; COI -Selected model: SYM+I+G. For reasons described above, only the model, not the parameters, has been fixed in the Bayesian analysis. Manual inspection of the likelihood scores of the 2630,001 trees samples in two parallel runs showed a good convergence after 650 trees, a value that has therefore been chosen as the burn-in. Posterior probabilities have been computed from a majority rule consensus tree of the remaining 2629,351 trees. The Bayesian tree (not shown, but values in Fig. 4), as well as the maximum likelihood tree (Fig. 4) show a topology identical to the maximum parsimony trees (Fig. 3). The branch length of the likelihood tree indicates larger intraspecific divergence in A. erythrocephalus, A. hova and A. jeekeli n. sp. than in all other analysed species (Fig. 4).

Biogeographical Patterns of Aphistogoniulus
Members of Aphistogoniulus are generally distributed in the rainforests and montane rainforests of Madagascar, spanning from the northern massif of Marojejy towards the south at Andohahela and Ambatotsirongorongo. Aphistogoniulus species are often tree-climbers [34] and show a widespread distribution (Fig. 1). The elevation-range, even on species-level, varies from close to sea level to mountaintops [34]. Only one Aphistogoniulus  species, A. aridus, is currently known from a dry habitat, the southern spiny forest of Mahavelo. In the far north of Madagascar, Aphistogoniulus seems to be replaced by its sistergenera Colossobolus and Sanguinobolus (Fig. 2). In this context, RASP indicates that the most likely ancestral habitat for the genus Aphistogoniulus was lowland rainforest (90.2%, Fig. 5) which is still inhabited by five species (Aphistogoniulus jeekeli, A. corallipes, A. infernalis, A. hova, and A. erythrocephalus). Interestingly, two species that are restricted to montane rainforests (Aphistogoniulus cowani and A. sanguineus) represent the sistergroup of two other species which are found in both lowland and mid-elevation rainforests (Aphistogoniulus erythrocephalus and A. hova). These results give evidence for a trend of colonization of higher rainforest levels from lowland rainforests within this clade. Two of the three Aphistogoniulus clades match a general distribution pattern; the A. erythrocephalus clade is distributed in the northern and median longitudes (Fig. 2, middle) in both lowland and mid-elevation rainforests. All representatives of the A. diabolicus clade can be found south of the 21u longitude (Fig. 2,  right), and include the only species of the genus recorded from an arid ecosystem. Inside the A. diabolicus clade, the northern-most species (A. jeekeli) is the most basal species (Fig. 3).
The distribution of the A. cowani clade is not consistent with this coherent biogeographical pattern. Its three members show a disjunctive distribution in mountains of the south (A. cowani), the middle (A. sanguineus) and the north (A. rubrodorsalis) of Madagascar (Fig. 2). The A. cowani clade also shows an unusual pattern regarding the altitude: its members seem to be restricted to montane forests, only living in areas of more than 900 m elevation.

Speciation Model for Aphistogoniulus
The majority of Aphistogoniulus species, including all members of the A. erythrocephalus and the A. diabolicus clades, show a stochastical distribution pattern which can best be explained by a 'simple' model reflecting speciation-by-distance, or simply coherent parapatric speciation. Our analyses give evidence that the most likely habitat for Aphistogoniulus was lowland rainforest (90.2%, Fig. 5), supporting these ideas. Major rivers do not appear to delineate species distribution borders (Fig. 1). At two nodes (Fig. 2), the phylogeny of Aphistogoniulus reflects two speciation models that are currently the subject of intense discussion [12].
(1) The A. cowani clade, clade D in our cladogram (Fig. 2), includes species restricted to areas with an elevation higher than 900 m. Here, the most congruent evolutionary model explaining the existence of such a clade is mountain refugia [13]. The ancestor of the A. cowani clade was probably a widespread species inhabiting lowland rainforests (63.6% probablity, see Fig. 5). The recent species are distributed in the Andringitra and Ranomafana massif (A. cowani), the Ambohijanahary and Ambohitantely Mountains (A. sanguineus), and Marojejy (A. rubrodorsalis). Distribution and speciation patterns reflecting the A. cowani clade were previously discovered in dwarf chameleons of the genus Brookesia [14], which occur in a similar ecological niche (leaf litter and low vegetation).
(2) The A. vampyrus -A. aridus species pair is an example of speciation possibly caused by an ecotone shift. The species A. vampyrus is like all other Aphistogoniulus species distributed in the rainforest, in this case the rainforests of Andohahela in southeastern Madagascar. Its sister-taxon A. aridus is currently only known from the Mahavelo forest, a spiny forest located approximately 60 km from Andohahela [51]. Here we have an ancestor (potentially A. vampyrus) originating in the south-eastern rainforest colonizing and adapting to the southern spiny forest (A. aridus). This speciation model was coined ''ecotone shift'' in other Malagasy animals, including amphibians and lemurs [12], while in Malagasy giant pill-millipedes such an ecotone shift seems to have occurred in the opposite direction from Aphistogoniulus: within the genera Zoosphaerium [52] and Sphaeromimus [53], the most basal species are distributed in the dry spiny forest ecosystem [10,11].

Molecular phylogeny of Aphistogoniulus vis-à-vis other molecular millipede studies
This is the first molecular millipede study on species level in the Southern hemisphere, and therefore the first on Madagascar. Millipede molecular studies on species level are still less than 10 and greatly focus on North America [35], with a single one also conducted in Japan [54]. In North America, Ice Age refugia were discovered to play a role in the current distribution of haplotypes of the spirobolid millipede genus Narceus. Unfortunately, the historical climate on Madagascar is still too little known for any comparison. There does not seem to be any indication of glaciation during the last ice ages on Madagascar. While the molecular Spirobolida dataset cannot be dated due to lack of fossils, the relatively large genetic distances between the different Aphistogoniulus species (14-30% of COI, 5-13% in the 16S, Supporting Information S5) are an indication that the speciation events in Aphistogoniulus predate the most recent ice ages.

Divergent Gonopod Evolution in the A. diabolicus Clade
The disc of the posterior gonopods is shaped like an almost closed 'C' in most members of the A. diabolicus clade, including A. jeekeli n. sp., A. diabolicus, A. vampyrus and A. aridus (Fig. 6), and displays only minor differences between species. However, the two other species of the A. diabolicus clade, A. corallipes and A. infernalis, show strongly modified posterior gonopods, where in the former the main branch and in the latter the basal branch are greatly enlarged (Fig. 6). This finding needs further investigation (especially if the corresponding female organs are also enlarged), and is another hint that divergent evolution in sexual morphology of related organisms with limited dispersal ability plays a role in speciation [54].

Analysis Problems
The current morphological dataset of Aphistogoniulus is comprised of only a very limited number of characters (Supporting Information S5). Although all higher clades are supported by good apomorphies, more characters would give a higher resolution to the analysis on the species level. While more morphological characters (such as the endotergum or the mandible) could be explored with a scanning electron microscopy study, it is questionable if (a) enough specimens of all rarely collected Aphistogoniulus species can be made available, and (b) if such characters really vary enough and are independent enough to provide a more robust phylogeny of Aphistogoniulus. Our reconstruction of the historical biogeography of Aphistogoniulus was also restricted to the molecular data set.
The three Aphistogoniulus species missing from the molecular analysis are those species only known from their type series. The absence of sequence data for A. rubrodorsalis n. sp., A. diabolicus and A. aridus is regrettable, but all three species are morphologically close enough to their respective sister taxa. Therefore, their exclusion from the molecular analysis does not likely affect the obtained results. Derivatio nominis. Rubrodorsalis, adjective, consisting of the assembled Latin words ruber = red and dorsalis = back, referring to the blood-red dorsal side of the species. Diagnosis. Main branch of telopodite of posterior gonopods with lateral membranous fringe strongly developed into erect process (Fig. 7B), a character only shared with A. sanguineus and A. cowani, identifying this species as a member of the A. cowani clade (Fig. 2, clade D). Distinguishable from the latter two species by a massive and apically curved basal branch of telopodite of posterior gonopods (Fig. 7D), which is slender in both A. cowani and A. sanguineus, and only in A. sanguineus more conspicuously curved. Retrorse process of anterior gonopod in A. rubrodorsalis n. sp. is more strongly developed than in A. sanguineus and A. cowani (Fig. 7C). Apical lobe of main branch of posterior gonopod folded backwards at a straight angle (Fig. 7D), while it is more diagonally folded in A. cowani and A. sanguineus.
Color (in preserved specimen). Head, antenna, legs and telson red. Anterior half of collum red, posterior half black. Body rings ventrally black, dorsally metazonites black and mesozonites red. Ozopores located in black area.
Conservation and distribution. A. rubrodorsalis is currently only known from the type locality, the Marojejy Mountain (Fig. 1 Diagnosis. Special shape and sharp-edged ending of basal branch of posterior gonopod (Fig. 8B), and presence of an apicomesal fringe on its main branch (Fig. 8D) identify A. jeekeli as a member of the A. diabolicus clade (Fig. 2, clade F). Membranous folds also present mesally on main branch of posterior gonopod (Fig. 8B), a character shared only by A. aridus and A. vampyrus. A. jeekeli can be distinguished from A. aridus and A. vampyrus by the presence of a unique rounded lobe (Fig. 7D) on the basal branch of posterior gonopod: only sclerotized spines are present in A. aridus and A. vampyrus. Main branch of posterior gonopod apically undulated (Fig. 8B). Coxite process of anterior gonopod slightly longer than telopodite (Fig. 8C). Collum completely red.
Color of head, antenna, collum, legs and telson red. Mesozonites of body rings dorsally red, ventrally darker. Anterior half of metazonite black, posterior half red. Red color extends ventrally along whole metazona. Ozopores marked in last third of body by black spots, which become larger towards telson.
Posterior gonopod telopodite branches forming an 'O' (Fig. 8B). Tips of main and basal branch close to but not touching one another. Main branch at mesal margin and apically with membranous folds, tip bi-lobed (Fig. 8D). Main branch slenderer and shorter than basal branch. Basal branch along its length of equal width, but tapering at apex (Figs. 8B, D); mesally with unique rectangular lobe (Fig. 8D).
Intraspecific variation. There is a difference in the number of body rings between the two observed populations. The specimens of Vatovavy possess 57 (F) and 58 (M) body rings, whereas the male specimens from Manombo only have 53 or 54 body rings and the two females 50 body rings, respectively. The population from Manombo is only tentatively placed in A. jeekeli because of strong morphological similarities, their gonopods and coloration are very similar. The large genetic distances between both populations (COI: 15.2%, 16S: 8.2%, 18S: 0.8%, see Supporting Information S10) are, however, only based on a single specimen from each population. Larger samples should be analyzed to discover if the population from Manombo might represent a distinct species.
Conservation and distribution. A. jeekeli is distributed in lowland rainforests of southeastern Madagascar. Both known localities, Vatovavy and Manombo, are isolated rainforest vestiges. Although the locations are 180 kilometers apart, all natural forests in between the locations are now degraded to pseudosteppe [54]. None of the current areas of distribution are effectively protected.

Outlook
The genetic differences between different populations of the widespread northern species A. erythrocephalus and A. hova are conspicuous enough to look into biogeographic patterns and possibly even cryptic speciation in those taxa when more properly preserved specimens become available. Further studies are necessary to highlight intraspecific differences or even cryptic speciation in two species of the A. cowani clade. Populations of A. cowani show color differences (amount of red on the collum, [34]) between populations in Ranomafana and Andringitra-Ivohibe (100 km distance). Also, the long distance (180 km) between both known populations of A. sanguineus (Ambohijanahary and Ambohitantely) needs to be further investigated for eventual gene flow and/or undiscovered populations in intermediate areas. Unfortunately, no appropriately conserved material from any of these populations is currently available for a genetic study.
The discovery of two new Aphistogoniulus species in the last year highlights how little we know about the millipede fauna of Madagascar. It is already difficult to say if the current distribution of the highland species A. sanguineus is so disjunctive because of recent habitat fragmentation or ecological specialization. Highland vegetation between the closely related A. sanguineus and the eastern species A. cowani and A. rubrodorsalis is now completely degraded [55,56], creating a void of information that hampers biogeographic studies [17].
Phylogenetic data remains lacking for most other millipede genera from Madagascar. In the order Spirobolida, only the genus Spiromimus of the subfamily Spiromiminae shares the wide distribution pattern of Aphistogoniulus [40], all other genera have a much more restricted range. A phylogenetic analysis of the Spiromiminae based on a morphological dataset showed several ecotone shifts occurring between rainforest and dry forest species of this subfamily [40]. The more mobile habit of the fire-millipede genus Aphistogoniulus, with its numerous widespread species, is in itself unusual for Malagasy millipedes. The wide range of the genus Aphistogoniulus and most of its species (Fig. 1), especially in comparison to other Spirobolida genera on Madagascar, could be related to its adaptation to a specific ecological niche, the surface of the leaf litter and lower branches of trees [34].
It might be interesting to compare the speciation pattern of the mobile Aphistogoniulus species with those of other Malagasy millipedes which occupy a different, more cryptic ecological niche, such as the wood-living genus Hylekobolus Wesener, 2009 (Spirobolida: Spirobolellidae [37]), or the giant pill-millipede genus Sphaeromimus (Sphaerotheriida: Arthrosphaeridae) that lives in the leaf litter.

Supporting Information
Supporting Information S1 New locality data and notes for described Aphistogoniulus species.