The structure of the male gonopods of millipedes has been considered to be species-specific. As such, gonopods—which aid in copulation and sperm transfer—are used in the taxonomic diagnosis and description of species. However, it was recently demonstrated that gonopod morphology is not always characteristic of species. Diagnoses based on gonopod morphology can therefore result in underestimation of taxonomic diversity amongst millipedes. On the basis of this observation, we examined genetic variation in two populations (approximately 250 km apart) of a widely distributed and colour-polymorphic southern African millipede, namely Bicoxidens flavicollis Attems, 1928. An analysis of genetic divergence based on 520 nucleotides of the mitochondrial cytochrome oxidase 1 gene, and 684 nucleotides of the mitochondrial 16S rRNA gene, demonstrated high levels of divergence (19.09% for cytochrome oxidase 1 and 6.66% for 16S rRNA) between the two populations. These results suggest the presence of cryptic species in B. flavicollis and, furthermore, corroborate observations that taxonomy based on gonopod morphology may be too inclusive.
Morphology — based taxonomy is extensively used to delimit taxa (Schlick-Steiner et al. 2007), particularly in groups that are poorly known and difficult to identify, such as millipedes (Hamer 2000). Male genitalia/gonopods are intromittent organs which are widely used in taxonomy on the grounds that they exhibit species-specific characters (Song & Bucheli 2010). Despite increasing use of DNA sequence data in invertebrate taxonomy (e.g. Pfenninger et al. 2007; Burns et al. 2008), gonopod morphology is still central to spirostreptid millipede taxonomy (e.g. Hoffman 2008; Hamer 2009; Mwabvu et al. 2010), because it has been mooted that the divergent male genitalia suggest reproductive isolation (Bond et al. 2003).
Although male gonopods evolve rapidly and divergently (Song & Bucheli 2010), speciation may be unaccompanied by change in gonopod morphology. According to Bond et al. (2003), speciation occurred without gonopod divergence in a species of the spirobolid millipede genus Anadenobolus Silvestri, 1897. In view of this evidence, morphology-based classifications are being re-evaluated against other criteria because morphology may fail to separate genetically distinct species (Bond & Sierwald 2002; Adams et al. 2009).
Sequences of mitochondrial cytochrome c oxidase 1 (CO1), 16S rRNAand 18S rRNA genes have been used to assess genetic divergence within many taxa. For example, mitochondrial 16S rRNA gene sequences were used in taxonomic studies of Chilopoda (Edgecombe & Giribet 2004), Polydesmida (Marek & Bond 2007), Hymenoptera (Dowton & Austin 1994), Australian elapid snakes (Keogh et al. 2000) and the Anadenobolus excisus (Karsch, 1881) millipede species complex (Bond & Sierwald 2002).
The genus Bicoxidens Attems, 1928, is endemic to southern Africa (Mwabvu et al. 2007). It includes medium to large (length 75–170 mm, diameter 6.6–10.6 mm) species (Mwabvu 2000; Mwabvu et al. 2007). Body colour ranges from shades of black and brown through to orange-yellow (Mwabvu et al. 2007). Bicoxidens species have been collected in diverse vegetation types, including Brachystegia or Acacia-dominated savannah woodland, and riverine and montane habitats (Mwabvu et al. 2007).
The taxonomic validity of Bicoxidens, as based on gonopod morphology, is not in doubt. Among the nine known species, B. flavicollis Attems, 1928, is the most widely distributed, occurring in east, north, south-west and central Zimbabwe, and in western Mozambique (Mwabvu et al. 2007) (Fig. 1). Although the male gonopods of B. flavicollis are identical in shape, body colour varies throughout the distribution range. Black, brown and orange-yellow specimens have been recorded in different habitats (Mwabvu et al. 2007). Given their poor dispersal ability and preference for moist microhabitats, millipedes are likely to speciate in isolation (Hopkin & Read 1992; Hamer & Slotow 2000). As a result, many millipede groups may contain hidden species (Brewer et al. 2012) that cannot be identified using gonopod morphology alone. Because reproductive isolation and genetic divergence are expected to increase with time, they are correlated (Fitzpatrick 2002). As such, the colour polymorphic and widely distributed B. flavicollis could be a species complex.
Considering the urgent need for information on invertebrate diversity in southern Africa (Hamer & Slotow 2002), there is much to be gained by investigating intraspecific genetic variation in millipedes, because diversity data and correct identification of taxa have implications for biodiversity conservation and other disciplines. Moreover, there is a paucity of genetic information on millipedes. This paper provides basic data and hopefully will stimulate further research on millipede genetics.
The objective of this research was to test the hypothesis that B. flavicollis is a species complex, by investigating genetic divergence in the CO1 and 16S rRNA mitochondrial genes in two populations. Given the great distance separating the two populations and the poor dispersal ability of millipedes, we predicted strong genetic differentiation between these populations of B. flavicollis.
MATERIAL AND METHODS
Fresh males of B. flavicollis were collected in Zimbabwe at Muterere (18°25′S 32°57′E) in the eastern highlands, and from Chihota, near Harare (18°15′S 31°05′E). In addition, males of Cacuminostreptus mazowensis Mwabvu in Mwabvu et al. 2010 obtained at Mazowe Dam (17°30′S 30°58′E) and Archispirostreptus tumuliporus (Karsch, 1881) from Marange (18°57′S 32°27′E), both localities being in Zimbabwe, were included as outgroups from the same family (Spirostreptidae). The millipedes were collected by hand, preserved in 100% ethanol, and deposited in the KwaZulu-Natal Museum, Pietermaritzburg, South Africa.
Total genomic DNA was obtained from the legs of preserved millipedes. The legs were first ground using a pestle and mortar. DNA was then extracted using the Qiagen DNeasy Blood and Tissue Kit (Qiagen) according to the manufacturer's instructions. Polymerase chain reactions targeting part of the mitochondrial CO1 and 16S rRNA genes were performed after optimisation using gradient PCR. Primers for the amplification of the genes were as used by Lavrov et al. (2002) and Bond & Sierwald (2002).
PCR products were electrophoresed at 15 V for 16 hours in 1.5% agarose gel containing 200 µl ethidium bromide (0.05 mg/m1) and using 0.5× TBE buffer. The PCR products were purified with a Zymoclean™ Gel Recovery Kit (Zymo Research, USA). Products were sequenced on an ABI 3730 capillary sequencer at Inqaba Biotechnical Industries (Pty) Ltd (Hatfield, Pretoria, South Africa). DNA samples that did not amplify after PCR, were cloned and sequenced for CO1. Raw sequences were edited in BioEdit version 7 (Hall 1999), and aligned using the Clustal W algorithm in BioEdit and by inspection. The alignments were trimmed to 684 nucleotides and 520 nucleotides of mitochondrial 16S rRNA and CO1 genes, respectively.
Bayesian, maximum parsimony and neighbour-joining analyses were performed to determine genetic distances between and phylogenetic relationships among the taxa. The general time-reversible (GTR) model was selected using the Akaike Information Criterion (AIC) for both the 16S rRNA and CO1 datasets, and was used subsequently in Bayesian and neighbour-joining analyses. Neighbour-joining and maximum parsimony analyses were implemented in PAUP (Swofford 2002). Genetic distances were presented as distance matrices and as a neighbour-joining tree, which was bootstrapped using 1000 pseudo-replicates. For parsimony analysis, the random additions sequence option (n=100) for discrete, unordered characters was used. The shortest tree was obtained using the heuristic search option under the tree bisection-reconnection (TBR) branch-swapping option. The degree of support for each node of the resulting tree was estimated using bootstrap re-sampling analysis (1000 pseudo-replicates; Felsenstein 1985). Bayesian analysis was implemented in Mr Bayes 3.0 (Huelsenbeck & Ronquist 2001), using flat priors. For all analyses, four Markov chains were run for 5 million generations each, and the first 500,000 trees were discarded as burn-in.
Comparisons of genetic divergence between Bicoxidens flavicollis specimens and the outgroups, indicated by general time-reversible genetic distances (%) based on 684 nucleotides of the mitochondrial 16S rRNAgene. Abbreviations: C — Chihota, D — Muterere, numbers next to each letter reflect the number of replicates.
Sequence data were registered at GenBank: KF057753 (B. flavicollis D5), KF057754 (B.flavicollis D3), KF057755 (B.flavicollis D1), KF057756 (B.flavicollis C5), KF057757 (B. flavicollis C2), KF057752 (C. mazowensis) and KF057758 (A. tumuliporus). Additional outgroup sequences of taxa belonging to the orders Spirobolida (Narceus americanus (Palisot de Beauvois, 1817), Narceus annularis (Rafinesque, 1820), Narceus gordanus (Chamberlin, 1943) and Anadenobolus excisus), Spirostreptida (Thyropygus sp.) and Julida (Cylindroiulus caeruleocinctus (Wood, 1864) and Anagaiulus blancatypa Enghoff, 1992) were obtained from the NCBI GenBank database.
We were able to obtain good quality sequence reads from only two specimens (one from each population) for CO1, but from three specimens from Muterere and two from Chihota populations for 16S RNA. The sequences from Chihota and Muterere were very divergent, and both CO1 and 16S genetic distances between the two populations were higher than expected. Genetic divergence between Chihota and Muterere specimens based on analysis of the CO1 gene, was 19.09 %. In respect of the 16S data, B. flavicollis was divided into two strongly supported sister clades (Fig. 2), separated by a genetic distance of greater than 6% (Table 1). Based on the 16S rRNA gene, the highest intrapopulation genetic variation at Chihota was 1.54 %. At Muterere, the range of variation was 0–0.15%, which was much lower than that at Chihota (Table 1). The 16S rRNA sequence dataset provided robust support for B. flavicollis from Chihota and Muterere as representing distinct taxa (Bayesian posterior probability 0.92 and 1.00, maximum parsimony bootstrap 100% and neighbour-joining bootstrap 100%) (Fig. 2).
Although only a small number of samples were successfully sequenced, genetic variation among specimens of B. flavicollis from the Muterere and Chihota populations nevertheless suggests that diagnosis based on gonopod morphology underestimates taxonomic diversity. These results lend support to the views of Hamer and Slotow (2000). After studying the gonopods of a common African millipede genus Doratogonus, Hamer and Slotow (2000) reported that the diagnoses of the species might be too inclusive. Studies in other invertebrate groups also support the assertion that morphology alone may not accurately reflect taxonomic diversity. For example, DNA sequences identified cryptic species in the lepidopteran genera Cymothoe Hüber (van Velzen et al. 2007) and Perichares Scudder (Burns et al. 2008), in an ant Tetramorium Mayr (Schlick-Steiner et al. 2006), and in the dipteran genus Chironomus Meigen (Pfenninger et al. 2007).
Although sequence divergence threshold values of 2–3 % have been assigned for designating species in insects and mammals (Rubinoff et al. 2006), the levels of sequence divergence and threshold values in African millipedes have not been determined. Therefore, this study is an important contribution to understanding millipede genetic diversity. The large CO1 and 16S genetic distances between B. flavicollis at Muterere and Chihota suggest the existence of more than one species. The contention of the presence of cryptic species in B. flavicollis would be supported by Lefébure et al. (2006), who proposed designation of species in Crustacea using a CO1 divergence threshold of below 16%, which is lower than the 19.09 % detected between the two populations of B. flavicollis. Moreover, with regard to 16S rRNA, genetic divergence of less than 2% is accepted as representing intraspecific variation in other taxa (Bond 2004). Therefore, if 2% divergence is used as a threshold value for species designation, 6% divergence between the two populations of B. flavicollis likewise suggests the presence of cryptic species, which could not be identified using genitalic morphology. The two strongly supported subclades of Bicoxidens (Fig. 2) further justify our assertion that the taxon B. flavicollis contains more than one species.
The marked genetic variability in B. flavicollis is not surprising. Bicoxidens species are small-bodied, with poor dispersal ability (Mwabvu et al. 2007), and as they tend to have very limited distributions, there is presumably restricted gene flow. Furthermore, the small population in Muterere occurs in vegetation on an isolated granite outcrop. This probably explains the low intra-population genetic variation compared with the population in Chihota. In addition, besides the genetic variation between the two populations, at Muterere the specimens of B. flavicollis are yellow with brown stripes, while at Chihota the millipedes are black, as in most areas of the distribution range.
Greater than expected genetic distinctiveness between the two populations of B. flavicollis emphasises the need to review taxonomies based on gonopod morphology. Because gonopod morphological change may not be correlated with genetic change (Adams et al. 2009), morphological ‘species’ might be hiding considerable levels of cryptic variation. Many species of millipedes are local or site endemics (Hamer & Slotow 2002), making them IUCN threatened (Hamer 2009), and thus in need of prioritising for conservation. Therefore, a threshold genetic distance for designating millipede species should be assigned in order to facilitate investigations into cryptic diversity.
In the light of differences in the rate of evolution of genes (Lefébure et al. 2006) and because of the small number of replicates in this study, the notion that B. flavicollis is a species complex warrants further research, using more genes and larger samples from throughout the distribution range. The information could be used to infer patterns and processes, and to enhance understanding of endemicity levels in spirostreptid millipedes in Africa. In addition to increased taxon sampling, further investigations into the suitability of the genes are required because some genetic markers are appropriate only for differentiating specific taxa (see Vences et al. 2005).
Given the large degree of genetic divergence between the two populations of B. flavicollis, the considerable intraspecific variation in gonopods reported in some groups (Brewer et al. 2012), and that reproductive isolation is correlated with genetic divergence (Fitzpatrick 2002), it would be valuable to investigate gonopod morphology in millipedes (as suggested by Brewer et al. (2012)) using morphometric landmark or shape analyses. Such studies might identify subtle variations in gonopods that demonstrate morphological stasis in shape and structure. Because the distribution of the nine currently recognised Bicoxidens species is known, future work could also include genetic divergence and cytogenetic studies of the species, in order to understand speciation in spirostreptids and to provide useful data on chromosomal evolution.
The research was funded by the University of KwaZulu-Natal and a Focus Area grant from the National Research Foundation, South Africa. We wish to thank Miss Theshnie Naidoo for training the first author in molecular laboratory techniques.
 http://africaninvertebrates.org urn:lsid:zoobank.org:pub:1BD03DEB-501E-4EF6-A0B8-A26A7AC99D36