Analyses of mitochondrial genes reveal two sympatric but genetically divergent lineages of Rhipicephalus appendiculatus in Kenya
Parasites & Vectors volume 9, Article number: 353 (2016)
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.
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.
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.
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 . 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 . For instance, different vertebrate hosts have been shown to influence the genetic structure of Ixodes uriae , 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 .
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 . 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 . 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 ). Heavy infestations can lead to tick worry, damaged hides - especially the ears where R. appendiculatus often congregate, anemia and toxicosis that results in enhanced susceptibility to other diseases .
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 , differences in body size [10, 14], vector competence  and in response to acaricides . 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–19]. Differences in agro-ecological and climatic conditions were thought to drive the differentiation of the two lineages [17–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 .
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 , 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 . 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.
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–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. . 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 ). COI gene was amplified using primers described in Folmer et al.  while the 12S rRNA gene was amplified using primers described in Simon et al. . 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  and reverse primer ITS2R . 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 MgCl2, 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 . 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 . 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–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”  following 1000 coalescent simulations. The mismatch distributions were augmented with the Fu’s F S  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 . The best nucleotide substitution model for the gene was T92 + G model  as determined with MEGA v6.0. Clade support 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  was constructed using COI sequences with NETWORK 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 haplotypes were deposited in the GenBank database under accession numbers KU725890-KU725917 while the protein identifiers for the corresponding translated protein sequences for these haplotypes were ANF89378-ANF89405. Two haplotypes (Hap_4 and Hap_1) were defined by 107 (32.2 %) and 94 (28.3 %) sequences, respectively, accounting for 60.5 % of all sequences analysed (Additional file 3: Table S3). Hap_4 was exclusive to Muguga low-line (ML), Lab OlPejeta (LP) and Zambia Eastern Province (ZE) populations and was observed in eight out of the 12 sequences from Kiambu High-line (KH). Hap_1 was exclusive to South Africa Lab (SAL), South Africa Natal (SAN), Kiambu unselected line (KU) and Zimbabwe West Mashonaland (ZM) and was also observed in eight of the 12 sequences from Muguga infected (MF) and Zambia Southern Province (ZS) and in 19 out of the 20 sequences from Kilifi (KF) (Table 2). All sequences from Uganda (UG) (n = 12) were exclusively of one haplotype (Hap_8). The average number of haplotypes across the 22 study populations was 28 and on average, 20, 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 , we reconstructed a ML tree using a 415 bp region derived from our 28 haplotypes combined with ten haplotypes defined by Mtambo et al. . We used a 415 bp fragment because this was the size of the fragment amplified by Mtambo et al. . 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 laboratory-bred 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 (Haplogroup A: SSD = 0.0126, P = 0.550; RI = 0.050, P = 0.740; haplogroup B: SSD = 0.0028, P = 0.170; RI = 0.102, P = 0.210). Both had negative Tajima’s D (-1.650 and -1.087, respectively) and Fu’s F S (-10.348 and−3.462, respectively) values (Table 4). While Tajima’s D (P = 0.017) and Fu’s F S (P <0.001) attained significance for haplogroup A, the values for haplogroup B did not (P = 0.125, P = 0.057) (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 . 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.
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–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–58]. It has also been used previously to show R. appendiculatus speciation [17, 18, 59] and the current study 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 . 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 .
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. , 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.  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 , physiological  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 , 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 sub-haplogroups 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.
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.
AMOVA, analysis of molecular variance; ARC-OVI, Agricultural Research centre- Onderstepoort Veterinary Institute; AusAID, Australian Agency for International Development; AWARD, African Women in Agricultural Research and Development; BecA, Biosciences eastern and central Africa; BLASTN, basic local alignment search tool; bp, base pair; CIDA, Canadian International Development Agency; COI, cytochrome c oxidase subunit I; CSIRO, Australia’s Commonwealth Scientific and Industrial Research Organisation; DAAD, German academic exchange service; DNA, deoxyribonucleic acid; dNTP, deoxyribonucletide triphosphate; dNTPs, deoxyribonucleotide triphosphates; ECF, east coast fever; ICARDA, International Centre for Agricultural Research in the Dry Areas; icipe, International Centre of Insect Physiology and Ecology; ILRI, International Livestock Research Institute; ITS2, Internal transcribed spacer 2; KWS, Kenya Wildlife Service; MJ, median joining; ML, maximum likelihood; mv, median vector; NCBI, National Center Bioinformatics; NEPAD, new partnership for Africa’s development; PCR, polymerase chain reaction; rDNA, ribosomal deoxyribonucleic acid; SFSA, Syngenta Foundation for Sustainable Agriculture; SSD, sum of squares deviation
Randolph SE. Ticks and tick-borne disease systems in space and from space. Adv Parasitol. 2000;47:217–43.
Price PW. Evolutionary biology of parasites. Princeton University Press; 1980
Poulin R, George-Nascimento M. The scaling of total parasite biomass with host body mass. Int J Parasitol. 2007;37:359–64.
Kain DE, Sperling FH, Daly HV, Lane RS. Mitochondrial DNA sequence variation in Ixodes pacificus (acari : ixodidae). Heredity. 1999;83:378–86.
McCoy KD, Boulinier T, Tirard C, Michalakis Y. Host-dependent genetic structure of parasite populations: differential dispersal of seabird tick host races. Evolution. 2003;57:288.
Cangi N, Horak IG, Apanaskevich DA, Matthee S, das Neves LCBG, Estrada-Peña A, et al. The influence of interspecific competition and host preference on the phylogeography of two African ixodid tick species. PLoS One. 2013;8:e76930.
Perry BD, Lessard P, Norval RAI, Kundert K, Kruska R. Climate, vegetation and the distribution of Rhipicephalus appendiculatus in Africa. Parasitol Today. 1990;6:100–4.
Randolph SE. Abiotic and biotic determinants of the seasonal dynamics of the tick Rhipicephalus appendiculatus in South Africa. Med Vet Entomol. 1997;11:25–37.
Walker JB, Keirans JE, Horak IG. The genus Rhipicephalus (Acari, Ixodidae): a guide to the brown ticks of the world. Cambridge: Cambridge University Press; 2005.
Norval RAI, Perry BD, Young AS. The epidemiology of theileriosis in Africa. London: London Academic Press; 1992.
Minjauw B, Mcleod A. Tick-borne diseases and poverty. The impact of ticks and tick-borne diseases on the livelihoods of small-scale and marginal livestock owners in India and eastern and southern Africa. UK: Vet. Med. DFID Animal Health Programme, Centre for Tropical Veterinary Medicine, University of Edinburgh; 2003.
Jongejan F, Uilenberg G. The global importance of ticks. Parasitology. 2005;129:S3–14.
Madder M, Speybroeck N, Brandt J, Berkvens D. Diapause induction in adults of three Rhipicephalus appendiculatus stocks. Exp Appl Acarol. 1999;23:961–8.
Shaw M, Young A. The biology of Theileria species in ixodid ticks in relation to parasite transmission. In: KF Harris editors. (Vol. 10) Adv Dis Vector Res New York, USA: Springer-Verlag New York Inc; 1994. p. 23–63.
Ochanda H, Young AS, Medley GF, Perry BD. Vector competence of 7 rhipicephalid tick stocks in transmitting 2 Theileria parva parasite stocks from Kenya and Zimbabwe. Parasitology. 1998;116(Pt 6):539–45.
Chigagure NN, Baxter GD, Stephen C. Microsatellite loci of the cattle tick Boophilus microplus (Acari : Ixodidae). Exp Appl Acarol. 2001;24(12): 951–6
Mtambo J, Madder M, Van Bortel W, Geysen D, Berkvens D, Backeljau T. Genetic variation in Rhipicephalus appendiculatus (Acari: Ixodidae) from Zambia: correlating genetic and ecological variation with Rhipicephalus appendiculatus from eastern and southern Africa. J Vector Ecol. 2007;32:168–75.
Mtambo J, Madder M, Van Bortel W, Chaka G, Berkvens D, Backeljau T. Further evidence for geographic differentiation in Rhipicephalus appendiculatus (Acari: Ixodidae) from Eastern and Southern provinces of Zambia. Exp Appl Acarol. 2007;41:129–38.
Speybroeck N, Madder M, Thulke HH, Mtambo J, Tirry L, Chaka G, et al. Variation in body size in the tick complex Rhipicephalus appendiculatus/Rhipicephalus zambeziensis. J Vector Ecol. 2004;29:347–54.
Madder M, Speybroeck N, Brandt J, Tirry L, Hodek I, Berkvens D. Geographic variation in diapause response of adult Rhipicephalus appendiculatus ticks. Exp Appl Acarol. 2002;27:209–21.
Kanduma EG, Mwacharo JM, Mwaura S, Njuguna JN, Nzuki I, Kinyanjui PW, et al. Multi-locus genotyping reveals absence of genetic structure in field populations of the brown ear tick (Rhipicephalus appendiculatus) in Kenya. Ticks Tick Borne Dis. 2016;7:26–35.
Radley DE. Infection and treatment method of immunization against theileriosis. In: Irvin AD, Cunningham MP, Young AS, editors. Advances in the control of theileriosis. Netherlands: Springer Publishing; 1981. p. 227–37.
Di Giulio G, Lynen G, Morzaria S, Oura C, Bishop R, Di Giulio G. Live immunization against east coast fever-current status. Trends Parasitol. 2009;25:85–92.
Young AS, Dolan TT, Mwakima FN, Ochanda H, Mwaura SN, Njihia GM, et al. Estimation of heritability of susceptibility to infection with Theileria parva in the tick Rhipicephalus appendiculatus. Parasitology. 1995;111:31–8.
Odongo DO, Ueti MW, Mwaura SN, Knowles DP, Bishop RP, Scoles GA. Quantification of Theileria parva in Rhipicephalus appendiculatus (Acari: Ixodidae) confirms differences in infection between selected tick strains. J Med Entomol. 2009;46:888–94.
Kanduma EG, Mwacharo JM, Sunter JD, Nzuki I, Mwaura S, Kinyanjui PW, et al. Micro- and minisatellite-expressed sequence tag (EST) markers discriminate between populations of Rhipicephalus appendiculatus. Ticks Tick Borne Dis. 2012;3:128–36.
Bailey KP. Noteson the rearing of Rhipicephalus aippendiculatus and their infection with Theileria parva for experimental transmission. Bull Epizoot Dis Afr. 1960;8:33–43.
Irvin AD, Brocklesby DW. Rearing and maintaining Rhipicephalus appendicutatus in the laboratory. Inst Anim Tech J. 1970;21:106-112
Hoogstraal H. African ixodoidea. VoI. I. Ticks of the Sudan (with special reference to Equatoria province and with preliminary reviews of the genera Boophilus, Margaropus and Hyalomma). 1956.
Walker J. Notes on the common tick species of East Africa. Nairobi: Copper, McDougall and Robertson Ltd; 1970.
Walker AR, Bouattour A, Camicas J-L, Estrada-Pena A, Horak IG, et al. Ticks of domestic animals in Africa: a guide to identification of species. Edinburgh (Scotland) Bioscience Reports. 2003;165-221.
Folmer O, Black M, Hoeh W, Lutz R, Vrijenhoek R. DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Mol Mar Biol Biotechnol. 1994;3:294–9.
Simons C, Frati F, Beckenbach A, Crespi B, Liu H, Flook P. Evolution, weighting, and phylogenetic utility of mitochondrial gene sequences and a compilation of conserved polymerase chain reaction primers. Ann. Entomol Soc Am 1994. p. 651–701.
Barker SC. Distinguishing species and populations of rhipicephaline ticks with its 2 ribosomal RNA. J Parasitol. 1998;84:887–92.
Domanico MJ, Phillips RB, Oakley TH. Phylogenetic analysis of Pacific salmon (genus Oncorhynchus) using nuclear and mitochondrial DNA sequences. Can J Fish Aquat Sci. 1997;54:1865–72. NRC Research Press Ottawa, Canada.
Librado P, Rozas J. DnaSP v5 : A software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009;25:1451–2.
Excoffier L, Lischer HEL. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour. 2010;10:564–7.
Rogers AR, Harpending H. Population growth makes waves in the distribution of pairwise genetic differences. Mol Biol Evol. 1992;9:552–69.
Harpending HC. Signature of ancient population growth in a low-resolution mitochondrial DNA mismatch distribution. Hum Biol. 1994;66:591–600.
Rogers AR. Genetic evidence for a Pleistocene population explosion. Evolution (N. Y). Soc Study Evol. 1995;49:608–16.
Fu YX. Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics. 1997;147:915–25.
Tajima F. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989;123:585–95.
Tajima F. The amount of DNA polymorphism maintained in a finite population when the neutral mutation rate varies among sites. Genetics. 1996;143:1457–65.
Tamura K, Stecher G, Peterson D, Filipski A, Kumar S. MEGA6: molecular evolutionary genetics analysis version 6.0. Mol Biol Evol. 2013;30:2725–9.
Tamura K. Estimation of the number of nucleotide substitutions when there are strong transition-transversion and G + C-content biases. Mol Biol Evol. 1992;9:678–87.
Bandelt HJ, Forster P, Rohl A. Median-joining networks for inferring intraspecific phylogenies. Mol Biol Evol. 1999;16:37–48.
Hebert PDN, Ratnasingham S, deWaard JR. Barcoding animal life: cytochrome c oxidase subunit 1 divergences among closely related species. Proc Biol Sci. 2003;270:S96–9.
Hebert PDN, Cywinska A, Ball SL, deWaard JR. Biological identifications through DNA barcodes. Proc Biol Sci. 2003;270:313–21.
Lohman DJ, Ingram KK, Prawiradilaga DM, Winker K, Sheldon FH, Moyle RG, et al. Cryptic genetic diversity in ‘widespread’ Southeast Asian bird species suggests that Philippine avian endemism is gravely underestimated. Biol Conserv Elsevier Ltd. 2010;143:1885–90.
Bucklin A, Steinke D, Blanco-Bercial L. DNA barcoding of marine metazoa. Ann Rev Mar Sci Annual Reviews. 2011;3:471–508.
Saitoh T, Sugita N, Someya S, Iwami Y, Kobayashi S, Kamigaichi H, et al. DNA barcoding reveals 24 distinct lineages as cryptic bird species candidates in and around the Japanese Archipelago. Mol Ecol Resour. 2015;15:177–86.
Čandek K, Kuntner M. DNA barcoding gap: reliable species identification over morphological and geographical scales. Mol Ecol Resour. 2015;15:268–77.
Tavares ES, Gonçalves P, Miyaki CY, Baker AJ. DNA barcode detects high genetic structure within neotropical bird species. PLoS One. 2011;6:e28543. Public Library of Science.
Clare EL, Lim BK, Fenton MB, Hebert PDN. Neotropical bats: estimating species diversity with DNA barcodes. PLoS One. 2011;6:e22648. Public Library of Science.
Song S, Shao R, Atwell R, Barker S, Vankan D. Phylogenetic and phylogeographic relationships in Ixodes holocyclus and Ixodes cornuatus (Acari: Ixodidae) inferred from COX1 and ITS2 sequences. Int J Parasitol. 2011;41:871–80.
Burger TD, Shao R, Barker SC. Phylogenetic analysis of the mitochondrial genomes and nuclear rRNA genes of ticks reveals a deep phylogenetic structure within the genus Haemaphysalis and further elucidates the polyphyly of the genus Amblyomma with respect to Amblyomma sphenodonti and A. Ticks Tick Borne Dis. 2013;4:265–74.
Burger TD, Shao R, Barker SC. Phylogenetic analysis of mitochondrial genome sequences indicates that the cattle tick, Rhipicephalus (Boophilus) microplus, contains a cryptic species. Mol Phylogenet Evol. 2014;76:241–53.
Dantas-Torres F, Latrofa M, Annoscia G, Giannelli A, Parisi A, Otranto D. Morphological and genetic diversity of Rhipicephalus sanguineus sensu lato from the New and Old worlds. Parasit Vectors. 2013;6:213.
Murrell A, Campbell NJ, Barker SC. Phylogenetic analyses of the rhipicephaline ticks indicate that the genus Rhipicephalus is paraphyletic. Mol Phylogenet Evol. 2000;16:1–7.
Gadagkar SR, Rosenberg MS, Kumar S. Inferring species phylogenies from multiple genes: concatenated sequence tree versus consensus gene tree. J Exp Zool B Mol Dev Evol. 2005;304:64–74.
Balmer O, Beadell JS, Gibson W, Caccone A. Phylogeography and taxonomy of Trypanosoma brucei. PLoS Negl Trop Dis. 2011;5:e961. Public Library of Science.
Murrell A, Campbell NJ, Barker SC. Recurrent gains and losses of large (84–109 bp) repeats in the rDNA internal transcribed spacer 2 (ITS2) of rhipicephaline ticks. Insect Mol Biol. 2001;10:587–96.
Norval R, Lightfoot C. Tick problems in wildlife in Zimbabwe. Factors influencing the occurrence and abundance of Rhipicephalus appendiculatus. Zimbabwe Vet J. 1982;13:11–20
Okello-Onen J, Tukahirwa EM, Perry BD, Rowlands GJ, Nagda SM, Musisi G, et al. Population dynamics of ticks on indigenous cattle in a pastoral dry to semi-arid rangeland zone of Uganda. Exp Appl Acarol. 1999;23:79–88.
We are grateful to the Kenya Wildlife Service (KWS) and Ol Pejeta conservancy (Kenya) management for granting permission to sample ticks from areas under their jurisdiction. We wish to thank Heloise Heyne (Agricultural Research Centre - Onderstepoort Veterinary Institute (ARC-OVI)) for assistance with tick identification and Milton Ochuka, International Livestock Research Institute (ILRI) for assistance with tick collection.
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.
RPB and RAS conceived, designed, coordinated and supervised the study. EGK carried out field work, molecular studies, data analysis and drafted the manuscript. JMM participated in data interpretation, analysis and results presentation. EK provided some field tick specimens while SM provided laboratory stocks and were both involved in morphological identification of samples. PWK was a research supervisor while NG and LMK provided technical assistance during the early stages of the project. RPB, RAS, JMM, EGK, PWK and NG read and reviewed the manuscript to its current form. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
Number of ITS2 and 12S rDNA sequences analysed from eight field and four laboratory populations of R. appendiculatus. (DOCX 14 kb)
PCR primer sequences for each gene marker amplified and their corresponding annealing temperatures. (DOCX 15 kb)
R. appendiculatus COI haplotypes showing variable sites and the number sequences in each haplotype. (DOCX 20 kb)
About this article
Cite this article
Kanduma, E.G., Mwacharo, J.M., Githaka, N.W. et al. Analyses of mitochondrial genes reveal two sympatric but genetically divergent lineages of Rhipicephalus appendiculatus in Kenya. Parasites Vectors 9, 353 (2016). https://doi.org/10.1186/s13071-016-1631-1