Analyses of mitochondrial genes reveal two sympatric but genetically divergent lineages of Rhipicephalus appendiculatus in Kenya

Background The ixodid tick Rhipicephalus appendiculatus transmits the apicomplexan protozoan parasite Theileria parva, which causes East coast fever (ECF), the most economically important cattle disease in eastern and southern Africa. Recent analysis of micro- and minisatellite markers showed an absence of geographical and host-associated genetic sub-structuring amongst field populations of R. appendiculatus in Kenya. To assess further the phylogenetic relationships between field and laboratory R. appendiculatus tick isolates, this study examined sequence variations at two mitochondrial genes, cytochrome c oxidase subunit I (COI) and 12S ribosomal RNA (rRNA), and the nuclear encoded ribosomal internal transcribed spacer 2 (ITS2) of the rRNA gene, respectively. Results The analysis of 332 COI sequences revealed 30 polymorphic sites, which defined 28 haplotypes that were separated into two distinct haplogroups (A and B). Inclusion of previously published haplotypes in our analysis revealed a high degree of phylogenetic complexity never reported before in haplogroup A. Neither haplogroup however, showed any clustering pattern related to either the geographical sampling location, the type of tick sampled (laboratory stocks vs field populations) or the mammalian host species. This finding was supported by the results obtained from the analysis of 12S rDNA sequences. Analysis of molecular variance (AMOVA) indicated that 90.8 % of the total genetic variation was explained by the two haplogroups, providing further support for their genetic divergence. These results were, however, not replicated by the nuclear transcribed ITS2 sequences likely because of recombination between the nuclear genomes maintaining a high level of genetic sequence conservation. Conclusions COI and 12S rDNA are better markers than ITS2 for studying intraspecific diversity. Based on these genes, two major genetic groups of R. appendiculatus that have gone through a demographic expansion exist in Kenya. The two groups show no phylogeographic structure or correlation with the type of host species from which the ticks were collected, nor to the evolutionary and breeding history of the species. The two lineages may have a wide geographic distribution range in eastern and southern Africa. The findings of this study may have implications for the spread and control of R. appendiculatus, and indirectly, on the transmission dynamics of ECF. Electronic supplementary material The online version of this article (doi:10.1186/s13071-016-1631-1) contains supplementary material, which is available to authorized users.


Background
Knowledge relating to the intra-and inter-population genetic structure and variability amongst parasitic populations is important in understanding the dispersal and transmission dynamics of the pathogens they transmit. Several factors, including climate, host diversity, degree of tolerance of host species and control and management practices affecting host behavior are all thought to influence spatial distribution patterns of ticks [1]. The interaction between ticks and their hosts could result in genetic adaptations and divergence that may ultimately lead to genetic differentiation and speciation in ticks. The host's physiological, behavioral and demographic variability may also influence the genetic landscape of ectoparasites with limited dispersal ability such as ticks [2,3]. Other factors that are thought to influence the genetic variability of ticks include host availability and migration, ecological requirements of juvenile and adult stages, and tick dispersal ability [4]. For instance, different vertebrate hosts have been shown to influence the genetic structure of Ixodes uriae [5], while the availability of suitable hosts to the juvenile stages of Hyalomma rufipes and Amblyomma hebraeum can influence the geographical distribution of the adult stages of these two ixodid ticks [6].
Rhipicephalus appendiculatus is a three-host tick species whose ability to survive in a particular locality is determined by climatic conditions [7,8] and it almost entirely depends on its hosts for dispersal. It is widely distributed in eastern, central and southern Africa [9,10]. It lays eggs off its hosts and uses more than one host at different life-cycle stages, specifically larval, nymphal and adult instars. Large numbers of both adult and immature ticks can be found on cattle, goats, African buffalo (Syncerus caffer), Waterbuck (Kobus ellipsiprymnus), Eland (Taurotragus oryx), Greater kudu (Tragelaphus strepsiceros) and other large bovids [9]. The larval and nymphal stages frequently infest lagomorphs e.g. the Cape hare (Lepus capensis). Rhipicephalus appendiculatus is of major economic importance as the vector of the protozoan parasite Theileria parva, which causes East coast fever (ECF) in cattle [11]. Rhipicephalus appendiculatus also transmits Theileria taurotragi to cattle from Eland (Taurotragus oryx) causing benign bovine theileriosis, Anaplasma marginale resulting in bovine anaplasmosis, the nairovirus inducing Nairobi sheep disease, and Rickettsia conorii resulting in tick typhus in humans [9]). Heavy infestations can lead to tick worry, damaged hidesespecially the ears where R. appendiculatus often congregate, anemia and toxicosis that results in enhanced susceptibility to other diseases [12].
Several studies suggest that phenotypic diversity exists between different populations of R. appendiculatus. These include diapause in R. appendiculatus in southern Africa, which has not been observed in east African populations [13], differences in body size [10,14], vector competence [15] and in response to acaricides [16]. Morphological, physiological, epidemiological and phylogenetic data have shown the existence of two groups of R. appendiculatus in southern and eastern Africa, which were thought to represent two phylogeographically differentiated lineages [13,[17][18][19]. Differences in agro-ecological and climatic conditions were thought to drive the differentiation of the two lineages [17][18][19][20]. A recent analysis of micro-and minisatellite markers showed an absence of geographic and host-associated genetic structuring amongst field populations of R. appendiculatus in Kenya [21].
Several populations of R. appendiculatus have been maintained as laboratory stocks for sporozoite production and as representatives of field genotypes. For example, the standard laboratory stock of R. appendiculatus (designated Muguga) has been used to produce the Muguga cocktail vaccine against T. parva [22,23]. Previously, analysis of the biology of laboratory stocks of R. appendiculatus revealed differences in infection rates [24], and susceptibility to -and efficiency of acquisition of T. parva [15,25]. Recent assessments using micro-and minisatellite markers revealed distinct genetic groups in laboratory stocks of R. appendiculatus which were less diverse than their field counterparts [21]. Selection, reproductive isolation and inbreeding were thought to have led to the differentiation in the laboratory stocks. However, this finding has not been investigated further using genetic markers targeting the mitochondrial genome.
While the distribution of R. appendiculatus in Africa is determined by ecoclimatic factors, the genetic variability within the species remains poorly investigated. To further assess the phylogenetic relationships between field and laboratory R. appendiculatus tick stocks, this study examined sequence variation at the cytochrome c oxidase subunit I (COI) gene, 12S rDNA and the bi-parentally inherited ribosomal nuclear ITS2 region. The phylogenetic relationships, demographic dynamics and the partition of genetic diversity and structure amongst populations of R. appendiculatus were investigated.

Tick samples
The study used tick samples that had previously been described in earlier studies on population genetics of R. appendiculatus [21,26]. Genomic DNA from a total of 332 individuals from ten field populations and 12 laboratory maintained stocks of R. appendiculatus were used to sequence the mitochondrial COI gene. From the 332 samples, a subset of 93 samples from 12 populations was used to sequence the 12S rRNA gene while 87 ticks from the same subset were used to sequence the nuclear ITS2 gene spacer (Additional file 1: Table S1). These samples were randomly selected to represent tick populations falling within the two major COI haplogroups observed in this study. Of the ten field populations, six (118 individuals) came from areas grazed exclusively by cattle, two (43 individuals) from areas grazed exclusively by wildlife, and another two (46 individuals) came from areas co-grazed by wildlife and cattle. A total of 125 individuals were sampled from 12 laboratory colonies, which had been bred and maintained as closed genetic stocks (see [27,28]). One laboratory stock was originally sampled in Uganda (n = 12), one in Zimbabwe (West Mashonaland; n = 12) and two in Zambia (Eastern Province; n = 12; Southern Province; n = 8); the remaining eight stocks were collected in Kenya. The ticks had been identified following standard morphological criteria [29][30][31]. Details of the area of origin of the ticks, population and sampling site characteristics and the population codes used are as previously described in Kanduma et al. [26]. A list of all the study populations is given in Table 1.

DNA extraction and PCR amplification
The DNeasy® Blood and Tissue Kit (Qiagen GmbH, Hilden, Germany) was used to extract genomic DNA following minor modifications to the protocol (see [26]). COI gene was amplified using primers described in Folmer et al. [32] while the 12S rRNA gene was amplified using primers described in Simon et al. [33]. The ITS2 region (1-1.25 kb) was PCR amplified as two fragments: a full-length fragment, plus an internal 721 bp fragment to ensure good sequence coverage. The full-length fragment was amplified with the forward primer 3SAF [34] and reverse primer ITS2R [35]. The sequences of the primers used to PCR amplify the COI, 12S rRNA and the nuclear ITS2 fragment and their corresponding annealing temperatures are shown in Additional file 2: Table S2. All PCRs were carried out in 50 μl volumes containing 1X PCR buffer (Promega), 0.125 μmol MgCl 2 , 0.1 μM of each dNTP, 0.25 pmol of each primer, 1.25 U of Taq DNA polymerase (Promega) and 50 ng of template DNA. The PCR cycling profiles involved an initial denaturation at 95°C for 5 min followed by 35 cycles of 94°C for 1 min, annealing for 1 min (see Additional file 2: Table S2 for annealing temperatures) and extension at 72°C for 90 s for COI and 2 min for 12S rDNA and ITS2, respectively. A final extension step at 72°C for 10 min completed the amplification. PCR products were purified using the QIAquick® PCR Purification Kit (Qiagen GmbH, Hilden, Germany) following the manufacturer's protocol. The products were sequenced directly using the BigDye Terminator v3.1 cycle sequencing chemistry on an ABI 3730 DNA Analyzer in accordance with the manufacturer's methods (Applied Biosystems, UK).

Sequence editing and multiple alignments
All sequence chromatograms were visually inspected and the sequences edited manually using the CLC Main Workbench 6.8.3 (CLC bio, Qiagen GmbH, Hilden, Germany). The sequences were then trimmed to remove low quality reads at the 5' and 3' ends. Consensus sequences for each gene were generated from the sequenced fragments. Prior to analyses, all sequences were trimmed to uniform sizes (COI, 558 bp; 12S rDNA, 345 bp; ITS2, 1149 bp). Multiple sequence alignments were performed for each gene using ClustalW2 in CLC Main Workbench. Species identity was investigated and confirmed via BLASTN searches on the NCBI database (http://blast.ncbi.nlm.nih.gov/Blast.cgi).

Genetic variation and structure
Sequences were collapsed into haplotypes, following multiple sequence alignments, using DnaSP v5.10.01 [36]. Genetic variation represented as nucleotide and haplotype diversity and mean number of nucleotide differences for the COI gene were calculated for each population, groups of populations and haplogroups using DnaSP. The partition of genetic variation within and among populations was assessed via nested analysis of molecular variance (AMOVA) using Arlequin v3.5 [37]. The groupings used in AMOVA were as follows: (i) one group composed of all sequences of R. appendiculatus; (ii) two groups of sequences, i.e. those from areas grazed exclusively by cattle vs those from areas co-grazed by cattle and wildlife; (iii) two groups of sequences, i.e. those from areas grazed exclusively by cattle vs those from areas grazed exclusively by wildlife; (iv) two groups of sequences, i.e. those from areas co-grazed by wildlife and cattle vs those from areas grazed exclusively by wildlife; (v) two groups of sequences, i.e. field stocks vs laboratory stocks; (vi) three groups of sequences defined on the basis of the host species, i.e. cattle vs mixed cattle-wildlife vs wildlife, respectively; and (vii) amongst the groups identified by the phylogenetic and median-joining network analysis.

Demographic dynamics and phylogenetic structure
Demographic dynamics were inferred from mismatch distribution patterns [38][39][40] of COI haplotypes as implemented in Arlequin. The goodness-of-fit of the observed pattern of mismatches from the one expected under neutrality was tested using the sum of squares deviation (SSD) and Harpending's raggedness index "RI" [39] following 1000 coalescent simulations. The mismatch distributions were augmented with the Fu's F S [41] and Tajima's D [42,43] statistics which are also coalescent-based estimators of selective neutrality. Their significance was tested with 1000 coalescent simulations in Arlequin.
Phylogenetic reconstruction was performed using the COI gene employing the Maximum Likelihood (ML) algorithm implemented in MEGA v6.0 [44]. The best nucleotide substitution model for the gene was T92 + G model [45] as determined with MEGA v6.0. Clade support South Africa Natal province (SAN) South Africa Lab stock (SAL) Zambia Sothern province (ZS)  [24] was assessed via 1000 bootstrap replicates. To provide further support for the ML analysis and reveal in greater detail, and therefore gain further insights into the phylogeny of R. appendiculatus, median-joining (MJ) network [46] was constructed using COI sequences with NET-WORK 4.6 software (fluxus-engineering.com).

Confirmation of the species identification
The 332 samples used in this study generated 558 bp of high quality consensus COI sequences. Their molecular identity was confirmed via BLASTN searches against the NCBI's non-redundant nucleotide sequence database. The BLASTN searches returned high values of sequence similarity (97-100 %) with those of archived R. appendiculatus (GenBank AF132833; KC503257 and DQ859261).

COI sequence diversity
The 558 bp fragment of COI revealed 30 polymorphic sites which defined 28 haplotypes (Additional file 3: Table S3).
All of the 28 Total 107 28 94 21 a Total number sequences from each of the populations that were analysed (−) indicates that no samples from that particular population were included Only haplotypes represented by more than 20 sequences are shown 14, 7 and 5 haplotypes were observed in tick populations sampled from areas grazed by cattle only, co-grazed by cattle and wildlife, grazed by wildlife only and in laboratory stocks, respectively ( Table 1). The highest number of haplotypes (ten) was observed in Kitale (KT) and Field OlPejeta (FP) while the lowest (one) was observed in nine laboratory stocks. Haplotype sequences of each of the 22 studied populations were deposited in GenBank under the Accession numbers KX276862-KX276944 (Table 1). The haplotype diversity ranged from 0.900 ± 0.161 (mean ± standard deviation) in Ruma (RUM2) to 0 in nine laboratory stocks with an average value of 0.802 ± 0.014 (Table 1). Amongst tick populations sampled from the areas grazed by different host species, those from areas grazed exclusively by wildlife had the highest haplotype diversity (mean 0.767 ± 0.0064) and the laboratory stocks had the lowest (mean 0.143 ± 0.029). The average nucleotide diversity was 0.0123 ± 0.00019 ranging from 0 in nine laboratory stocks to 0.010 ± 0.06 in Kitale (KT) ( Table 1). The average number of nucleotide differences was 6.865 ± 3.2391 and ranged from 0 in nine laboratory stocks to 8.100 ± 4.534 in RM (Table 1). In general, ticks from field populations showed the highest levels of diversity whereas the laboratory stocks were the least diverse.

Phylogenetic relationships and median-joining network of COI haplotypes
To gain insights into the phylogenetic relationships between the 28 COI haplotypes, a ML tree (Fig. 1) and a MJ network were constructed (Fig. 2). The ML tree revealed two well-resolved groups of R. appendiculatus (bootstrap value of 100 %). The MJ network also revealed two groups that were separated by 12 mutation steps. The cluster of haplotypes in the ML tree and MJ network did not differ between the two algorithms. We therefore designated the two groups as haplogroups A and B, respectively. Haplogroup A clustered 19 haplotypes including Hap_4, the haplotype with the highest frequency, whereas haplogroup B contained nine haplotypes, which included Hap_1, the haplotype with the second highest frequency. Two median vectors (mv) were observed among the two haplogroups (Fig. 2); they may represent either haplotypes that were not sampled, or alternatively never present in Kenya, or have become extinct. A star-like pattern anchored by haplotypes H_4 and H_1 was evident for haplogroup A and B (Fig. 2), respectively, hinting at population expansion from an ancestral group, although the timescale is unclear. The ML tree (Fig. 1) appears to suggest the presence of two sub-haplogroups within haplogroup A (bootstrap value of 98 %). These can also be observed within the MJ network but are separated by a single mutation step. This suggests the possibility of genetic divergence within haplogroup A, requiring further analysis using a larger set of samples.
To test if the COI haplotypes generated in our study clustered with those of R. appendiculatus populations from eastern and southern Africa, which were also separated into two distinct groups [17], we reconstructed turanicus (JQ737086) from the GenBank database and another from a Kenya tick confirmed to be Rhipicephalus evertsi were included as the outgroup a ML tree using a 415 bp region derived from our 28 haplotypes combined with ten haplotypes defined by Mtambo et al. [17]. We used a 415 bp fragment because this was the size of the fragment amplified by Mtambo et al. [17]. Our haplotypes of haplogroup A clustered together with representative haplotypes from Zambia's eastern province and Rwanda whereas those of haplogroup B clustered together with representative haplotypes from the Comoro Islands and one haplotype each from Zambia's southern and eastern provinces, respectively (Fig. 3). Further examination of this tree reveals that the haplotypes that formed haplogroup A were subdivided into three sub-haplogroups (bootstrap values > 89 %) (Fig. 3). One sub-haplogroup (sub-haplogroup II) contained Kenyan haplotypes only (n = 5), another (sub-haplogroup I-B) comprised nine Kenyan haplotypes and one from Rwanda, while the third (sub-haplogroup I-A) was made up of five haplotypes from Kenya, two from Rwanda and four from Zambia's eastern province. This result suggests higher variation in R. appendiculatus, especially in haplogroup A, and a higher degree of phylogenetic complexity in this haplogroup not revealed in the studies of Mtambo et al. [17,18].

Population structure and demographic dynamics deduced from COI sequences
The partition of genetic variation within and among populations and groups of populations was investigated using hierarchical AMOVA taking into account seven population groups defined a priori ( Table 3). The highest level of genetic variation (90.8 %) was attributable to genetic differences between haplogroups A and B. Only 4.92 % of the total variation could be assigned to differences associated with the three host species complexes (cattle, cattle-wildlife and wildlife, respectively). Generally, the lowest levels of genetic variation were observed between groups of populations and ranged from−1.43 to 14.91 %. The variation present among individuals within populations ranged between 9.3 and 52.37 % whereas that among populations within groups was greater than 35 % in four comparisons. The observed variation between individuals within populations was greater than 43 % with the exception of the comparison amongst haplogroups (Table 3).
To provide insight into demographic dynamics, we analysed mismatch distribution patterns for different groups of populations. The overall mismatch distribution pattern for the 22 populations (Fig. 4a) was bimodal. The observed pattern did not deviate significantly from that expected under a model of expansion (SSD = 0.076, P = 0.09) and had a smooth distribution (RI = 0.061, P = 0.070) ( Table 1). The Tajima's D statistic was positive while Fu's F S was negative and neither were significant (Table 1). Taken together, these data are consistent with population expansion. We also investigated the demographic profiles for the ten populations sampled from the field and the 12 laboratory stocks (Fig. 4b, c). Both groups of populations exhibited two peaks. The observed pattern for the field populations deviated significantly from the one expected under a model of expansion (SSD = 0.845, P < 0.0001) with no significant variation around the curve (RI = 0.0317, P = 1.000). For this group, Tajima's D statistic was positive (D = 0.767) but this was not significant (P = 0.819), while Fu's F S parameter was negative (F S = -0.959) and not significant (P = 0.447). For the laboratory stocks, the observed pattern also deviated significantly from that expected (SSD = 0.189, P = 0.04). For this group of ticks, both Tajima's D statistic (D = 3.199, P = 0.996) and Fu's F S parameter (F S = 15.042, P = 0.994) were positive but not significant. This suggests that either the field populations have a weak signal of expansion or are in demographic equilibrium, whereas the laboratorybred stocks have been subject to an anthropogenic bottleneck and/or genetic drift. The two peaks observed in the overall dataset and in the field and laboratory stocks, respectively suggested the existence of two groups of ticks. The two peaks were found to correspond to the two haplogroups revealed by the ML and MJ analyses. We therefore performed mismatch analysis for each haplogroup (Fig. 4d, e). Both exhibited a unimodal profile and the observed patterns did not deviate significantly from that expected under a scenario of population expansion  Table 4). This suggests a strong signal of expansion for haplogroup A and a weaker one for haplogroup B. Our recent study utilizing nuclear satellite markers had also observed population expansion in field ticks [21]. These findings together with the star-like pattern observed in the MJ network, the mismatch distribution patterns and the two coalescent-based estimators of neutrality indicate expansion in the two haplogroups even in the absence of molecular dating.

Diversity and phylogenetic relationships based on 12S rRNA and ITS2 region
The 12S rRNA gene and ITS2 region were sequenced from a subset of the 332 R. appendiculatus individuals sequenced for the COI gene. Of the 93 12S rDNA sequences from 12 populations, five haplotypes were observed, two main (one defined by 38 sequences and the other by 52 sequences, respectively), and three minor (each defined by one sequence). Following ML phylogeny analysis, the five haplotypes clustered into two haplogroups which were identical to those generated from the COI gene.
A 1149 bp fragment of the ITS2 region was amplified from 87 individuals derived from different mitochondrial haplotypes. Three haplotypes were observed. One contained 67 sequences and the other two contained nine and 11 sequences, respectively. These ITS sequences did not cluster into groups corresponding to the COI or 12S rDNA haplogroups. The five 12S rDNA haplotype sequences were deposited in the GenBank database under accession numbers KX276945-49 and those of the three ITS2 haplotypes unde accession numbers KX276950-52.

Discussion
This study assessed the genetic relationships between populations of R. appendiculatus found in Kenya through the analysis of the mitochondrial COI and 12S rRNA genes and the nuclear transcribed ribosomal ITS2 fragment. COI gene has been and continues to be widely used as a marker for DNA barcoding to discriminate between closely related taxa [47][48][49][50][51][52]. Evolution of the COI gene is thought to be rapid enough to allow the discrimination of closely related species, as well as to detect intraspecific differentiation of phylogeographically distinct groups [53,54]. The utility of COI as a phylogenetic marker for ticks has been demonstrated previously [55][56][57][58]. It has also been used previously to show R. appendiculatus speciation [17,18,59] and the current study Clusters were based on a priori groupings of sampling localities. Cattle R. appendiculatus populations were collected directly from cattle or pastures grazed by cattle only. Cattle vs. cattle-wildlife and wildlife refers to populations collected from areas grazed by cattle versus a combination of populations from pastures co-grazed by cattle and wildlife and areas grazed by wildlife. Cattle vs wildlife only R. appendiculatus populations refer to ticks collected from areas grazed by cattle versus those collected areas grazed by wildlife only. Cattle-wildlife vs wildlife populations refer to populations from areas co-grazed by both cattle and wildlife versus wildlife only populations. Field vs laboratory R. appendiculatus populations refers to all R. appendiculatus ticks collected from field localities versus laboratory R. appendiculatus. Haplogroup A vs haplogroup B was between the two major R. appendiculatus haplogroups identified by ML and MJ network  found the variation in the COI to be adequate for phylogeny reconstruction and associated analyses using R. appendiculatus samples from Kenya.
In reconstructing the phylogenetic history of a species, the use of multiple genetic markers targeting different regions of the genome, is a better strategy in order to overcome the drawbacks of using a single marker, while increasing the accuracy of inference [60,61]. Here in addition to the COI gene, we analysed the phylogenetic relationships using the mitochondrially-encoded 12S rRNA gene and the nuclear genome-encoded ITS2 fragment. The COI analysis identified 28 haplotypes in 332 sequences. The NJ and MJ network partitioned these haplotypes into two distinct haplogroups. These two haplogroups were also discriminated by the 12S rDNA sequences but not by the nuclear transcribed ITS2 sequences. Using COI and 12S rDNA, Mtambo et al. [17,18] also observed two haplogroups of R. appendiculatus in eastern and southern Zambia but these were not detected by the ITS2 sequences. The low resolution afforded by ITS2 has also been reported in Amblyomma hebraeum and Hyalomma rufipes [6]. These findings suggest that COI and 12S rRNA genes are better markers for studying intraspecific diversity whereas the ITS2 fragment may be more useful in discriminating between species because it tends to show little intraspecific, but, considerable interspecific variation, possibly due to sexual recombination within species [62].
From the analysis of 332 COI sequences of R. appendiculatus, the overall mean number of nucleotide differences was 6.8647 ± 3.2391 and the mean haplotype and nucleotide diversities were 0.802 ± 0.014 and 0.0123 ± 0.0064, respectively. Cangi et al. [6], observed a lower level of haplotype and nucleotide diversities of 0.66 and 0.002, respectively, in A. hebraeum, an ixodid tick with a wider vertebrate host range, but a comparable level of haplotype and nucleotide diversity among isolates of 0.96 and 0.009, respectively, relative to the much more host-specialized H. rufipes. We also observed a high level of intra-and inter-population genetic diversity among the study populations. The values were much higher in the field ticks compared to the laboratory stocks, which were, by definition, subject to founder effects and population bottlenecks. The high diversity in field ticks is most probably the result of admixture between different geographic populations facilitated by the translocation of domestic animals either as trade items or through exchange following socio-cultural traditions. Indeed, no phylogeographic structure was revealed between the R. appendiculatus populations analysed in this study as revealed by either ML or MJ network analysis. In an earlier study, Kanduma et al. [21] observed no phylogeographic structure in field ticks that were analysed using autosomal micro-and minisatellite markers. The results suggest extensive translocation of ticks over a wide geographic range, in spite of low intrinsic dispersal ability of these arthropods resulting in populations with admixed genotypes. Domestic cattle in Kenya are frequently moved over large distances for commercial and socio-cultural reasons, as well as for seeking pasture during dry seasons. These would facilitate tick dispersal over a large geographical range, while the movement of the natural reservoirs of R. appendiculatus (wild bovidae) within the wildlife areas considered in this study is limited since these areas are fenced.
The laboratory stocks investigated here have been maintained as closed populations for over 30 years. It is therefore not surprising that they exhibited low levels of genetic diversity due to inevitably high levels of inbreeding. In spite of their inbred status, AMOVA revealed a negative value of genetic differentiation between the field and laboratory stocks implying that the two groups are much more related than might be expected. There are several potential explanations. First, that the inbreeding in the laboratory stocks has not resulted in a drastic reduction in their allelic variation; secondly, that variation present in the laboratory stocks is well represented in the field stocks; and thirdly, the induced bottleneck and genetic drift which could be due to inbreeding and small effective population sizes have not altered drastically their allelic composition.
Morphological [19], physiological [13] and phylogenetic [17,18] data previously identified two distinct groups of field R. appendiculatus in some parts of Africa and it was suggested that they may represent geographically differentiated lineages, that may have diversified as a result of distinct selective pressures. For instance, ticks found in southern Africa (South Africa, southern Zambia and Zimbabwe) and those found in eastern Africa (Kenya, Tanzania, Uganda, Burundi and Rwanda) were thought to constitute two geographically isolated groups of ticks that can be discriminated based on morphological, ecological and epidemiological differences [17,18]. In the current study, we observed two major haplogroups of R. appendiculatus in Kenya as defined by mitochondrial haplotype. These two haplogroups however exhibited no phylogeographic structure or correlation with the type of host species from which the ticks were collected or the evolutionary and breeding history of the species (field populations relative to laboratory stocks). Although we did not estimate the divergence time between the two genetic groups, it is possible that their divergence is not recent because they were observed among inbred laboratory stocks which were initially collected from populations of field ticks up to 50 years ago. AMOVA revealed that 90.8 % of the total genetic variation was explained by divergence within the two major haplogroups. While different host species have been shown to influence the spatio-genetic structure of other tick species, such as Ixodes uriae [5], the genetic variation between R. appendiculatus collected from different mammalian hosts was low (4.94 %). By contrast, the between population variation exceeded 35 %, whereas the variation between individuals within populations ranged between 9.3 and 52.3 %. This demonstrates low genetic differentiation between populations of R. appendiculatus sampled from different hosts suggesting minimal host specialisation. This suggests that genetic differentiation amongst tick populations in Kenya is a phenomenon primarily of ancestral differentiation between the two haplogroups and that recent reproductive isolation and the exploitation of different mammalian hosts has, to date, played a relatively minor role in driving this differentiation. The fact that the two major haplogroups that we have identified clusters together with representative haplotypes of R. appendiculatus from southern Africa [17,18], suggests a wide geographic distribution range of these haplogroups in eastern and southern Africa. It is possible that the original divergence in this species could have arisen either due to genetic drift and/or novel adaptations via selection giving rise to significant morphological, physiological and phenotypic changes seen in ticks from different geographical areas. Whether there are any associated phenotypic differences that can be used to discriminate the two haplogroups, which might influence parameters such as T. parva transmission dynamics, requires further investigation. Further investigation is also required into origin and possible evolutionary forces driving the occurrence of multiple subhaplogroups within haplogroup A.
We investigated the demographic dynamics of R. appendiculatus in Kenya by assessing the mismatch distribution patterns for the overall dataset, the field and laboratory stocks and within the two haplogroups identified by the ML and MJ network analysis. The results for the overall dataset, field populations and the two haplogroups suggest that these three groups of R. appendiculatus have passed through a demographic expansion perhaps associated with range expansion of a founder population.
The findings of this study may have taxonomic implications and suggest the potential for incipient speciation in R. appendiculatus. Rhipicephalus appendiculatus is a generalist tick, although Cape buffalo is the main wild host reservoir, and cattle are the preferred domestic hosts of the adult and nymphal instars [63,64]. Such a generalist ectoparasite which infests other wild and domestic animals can disperse across ecosystems potentially modifying disease transmission cycles. In this respect, understanding the population structure of R. appendiculatus is important in the design of sustainable control strategies, since different tick populations may be characterised by differences in vector competence, acaricide resistance and susceptibility to infection with T. parva. In future it will be important to establish how the phenotypes of the two R. appendiculatus haplogoups identified in this study differ, particularly with respect to acquisition and transmission of ECF.

Conclusions
COI and 12S genes are superior genetic markers for intra-species population genetic studies in R. appendiculatus over the ITS2. Based on these two genes, two distinct and well-differentiated haplogroups which have passed through a demographic expansion perhaps associated with range expansion of a founder population exist in Kenya. These two haplogroups have no phylogeographic structure or correlation with their mammalian host species or the evolutionary and breeding history of the species. There is a wide geographical distribution range of these two haplogroups in eastern and southern Africa. These findings may have important taxonomic implications and may point to an ongoing speciation of R. appendiculatus in sub-Saharan Africa. It would be important to establish if the two haplogroups have any associated phenotypic differences which might influence parameters such as T. parva acquisition and transmission dynamics. In addition, identifying evolutionary forces driving the observed genetic differentiation may help explain the apparent population expansion of the two haplogroups within the sub-Saharan region.

Additional files
Additional file 1:

Funding
This work was financially supported by the Biosciences eastern and central Africa Network (BecANet) through the New Partnership for Africa's Development (NEPAD), which was funded by the Canadian International Development Agency (CIDA) through a PhD fellowship to EGK. The study was also partially supported by the German Academic Exchange Service (DAAD) through a University of Nairobi PhD research grant to EGK. The African Women in Agricultural Research and Development (AWARD) also supported data analysis, thesis write-up and conference participation through a fellowship to EGK. We also gratefully acknowledge the financial support provided to the Biosciences eastern and central Africa Hub at the International Livestock Research Institute (BecA-ILRI Hub) by the Australian Agency for International Development (AusAID) through a partnership between Australia's Commonwealth Scientific and Industrial Research Organisation (CSIRO) and the BecA-ILRI Hub; and by the Syngenta Foundation for Sustainable Agriculture (SFSA), which made data interpretation, analysis and the thesis write-up possible.

Availability of data and materials
The dataset(s) supporting the conclusions of this article are available in the GenBank repository (http://www.ncbi.nlm.nih.gov/genbank/). The 28 haplotype sequences are under the accession numbers KU725890-KU725917 while their protein identifiers for the corresponding translated protein sequences have the numbers ANF89378-ANF89405. Haplotype sequences of each of the 22 studied populations are under the accession numbers KX276862-KX276944. The GneBank accession numbers of the five 12S rDNA haplotype sequences are KX276945-49 and those of the three ITS2 haplotypes sequences are KX276950-52.