Temporal stability of Glossina fuscipes fuscipes populations in Uganda.

BACKGROUND: Glossina fuscipes, a riverine species of tsetse, is the major vector of human African trypanosomiasis (HAT) in sub-Saharan Africa. Understanding the population dynamics, and specifically the temporal stability, of G. fuscipes will be important for informing vector control activities. We evaluated genetic changes over time in seven populations of the subspecies G. f. fuscipes distributed across southeastern Uganda, including a zone of contact between two historically isolated lineages. A total of 667 tsetse flies were genotyped at 16 microsatellite loci and at one mitochondrial locus.
RESULTS: Results of an AMOVA indicated that time of sampling did not explain a significant proportion of the variance in allele frequencies observed across all samples. Estimates of differentiation between samples from a single population ranged from approximately 0 to 0.019, using Jost's DEST. Effective population size estimates using momentum-based and likelihood methods were generally large. We observed significant change in mitochondrial haplotype frequencies in just one population, located along the zone of contact. The change in haplotypes was not accompanied by changes in microsatellite frequencies, raising the possibility of asymmetric mating compatibility in this zone.
CONCLUSION: Our results suggest that populations of G. f. fuscipes were stable over the 8-12 generations studied. Future studies should aim to reconcile these data with observed seasonal fluctuations in the apparent density of tsetse.

Introduction

Tsetse flies, Glossina spp (Diptera: Glossinidae) transmit several species of pathogenic trypanosomes causing Human African Trypanosomiasis (HAT) and African Animal Trypanosomiasis (AAT). HAT affects human welfare directly through the chronic and acute forms of the disease caused by Trypanosoma brucei gambiense and T. b. rhodesiense respectively. AAT, on the other hand, stands as a major obstacle to the development of more efficient and sustainable livestock production systems in tsetse-infested areas [1]. A major challenge to controlling HAT is lack of suitable prophylactic drugs and vaccines against trypanosomiasis. Furthermore, chemotherapeutic agents for treatment of HAT are expensive, difficult to administer in remote areas and exhibit poor safety profiles. Consequently, vector control remains a viable alternative for large-scale control of trypanosomiasis.

Understanding tsetse population dynamics is critical for determining which control strategy is most appropriate (e.g., suppression, eradication), for choosing the best method for enacting that strategy (e.g., traps, insecticide-treated cattle, sterile insect technique), and for determining the scale at which vector control activities must be implemented [2]. Determinants of population dynamics include both life history and ecological correlates such as mating system, dispersal ability and population size, which influence the extent to which tsetse populations can recover from refugia following intervention, or re-colonize a cleared zone from neighboring sources. Recently, the use of population genetics has provided insights into tsetse ecology [3], with important ramifications for the implementation of control programs [4]. For example, studies of tsetse in Guinea and Senegal have identified populations that are sufficiently isolated to warrant attempts at complete elimination [5-7]. Elsewhere though, studies have documented relatively high levels of gene flow, necessitating integration of barriers into elimination schemes [8] or warranting an area-wide control effort that encompasses the dispersal-linked populations [9,10].

Across Africa, Glossina fuscipes is one of the most important vectors of HAT, transmitting an estimated 90% of all disease cases [11]. Glossina fuscipes is a member of the palpalis group of tsetse, which inhabit low bushes or forests at the margins of rivers, lakes or temporarily-flooded scrub land. In eastern Africa, populations of the subspecies G. f. fuscipes appear to respond to seasonal weather patterns, often disappearing during the bi-annual dry season from sites where they were previously abundant [12]. If populations in refugia are small, then seasonal bottlenecks could result in large temporal changes in gene frequencies. In order to investigate the impact of seasonal climate changes on population size and to gain further insight into the population dynamics of G. f. fuscipes, we evaluated temporal changes in gene frequencies at one mitochondrial locus and 16 microsatellite loci in multiple Ugandan populations. Our sampling scheme included three populations situated at a zone of contact between two divergent lineages of G. f. fuscipes. These two lineages exhibit distinct mitochondrial DNA (mtDNA) haplotypes and strong differentiation at microsatellite loci, suggesting a long history of isolation, and providing a unique opportunity to monitor their interaction over time [9,10].

Materials and methods

Tsetse collection and study area

Tsetse flies were collected using biconical traps [13] during the period from March 2008 to January 2010. All sites were sampled in 2008 [10] and then at least one year later in 2009 or 2010. Four sites were also sampled a third time (Table 1). Each fly was stored individually in 80% ethanol.

Table 1

Indices of molecular diversity at mitochondrial and microsatellite loci for temporal samples of G. f. fuscipes.

		Microsatellites				Mitochondrial DNA
Sample	Date of Sampling	N	Allelic Richness	Ho	He	N	No. haplotypes	Haplotype diversity	Nucleotide diversity
BN - 0	March 2008	32	4.2	0.529	0.578	15	3	0.648	0.00538
BN - 8	March 2009	40	3.9	0.568	0.609
BN - 12	October 2009	64	4.1	0.549	0.574	16	4	0.692	0.00466
BU - 0	March 2008	39	3.5	0.459	0.485	17	1	0.000	0.00000
BU - 8	March 2009	40	3.4	0.476	0.485
BU - 12	October 2009	40	3.4	0.464	0.477	19	1	0.000	0.00000
JN - 0	March 2008	40	3.2	0.479	0.489	19	3	0.444	0.00731
JN - 13	January 2010	18	3.1	0.460	0.485	18	1	0.000	0.00000
MK - 0	March 2008	40	2.9	0.487	0.460	21	2	0.495	0.00093
MK - 8	March 2009	24	3.0	0.455	0.431
MK - 12	November 2009	22	3.1	0.418	0.445	21	2	0.467	0.00088
MS - 0	March 2008	40	3.7	0.568	0.547	18	2	0.471	0.00886
MS - 13	January 2010	17	4.4	0.562	0.597	17	3	0.559	0.00964
OK - 0	March 2008	39	3.3	0.452	0.507	17	3	0.471	0.00094
OK - 8	March 2009	40	3.4	0.563	0.546
OK - 12	October 2009	39	3.4	0.547	0.552	18	2	0.294	0.00055
OT - 0	July 2008	53	4.0	0.508	0.535	20	3	0.426	0.00122
OT - 11	November 2009	40	3.7	0.514	0.540	20	4	0.537	0.00131

Tsetse collections were conducted at seven sites spanning central and southeastern Uganda (Figure 1). These sites generally reflected the riverine/woodland habitat preferred by G. fuscipes, but varied somewhat in regard to the immediate environment. Sites at Busime (BU) and Junda (JN) were located in a transition zone between marsh and woodland on the edge of Lake Victoria and Lake Kyoga, respectively. Sites at Bunghazi (BN), Dokolo/Otuboi (OT) and Okame (OK) were situated along permanent streams in a region of mixed agriculture and pastureland. Sampling at Mukongoro (MK) was conducted at the margin of ephemeral wetlands associated with rice cultivation. Sampling at Masindi (MS) was conducted within a region of banana and sugar cane plantations.

Figure 1

Map of sites at which populations of . Location codes are shown in Table 1. The dotted line indicates the approximate extent of a zone of contact between two historically isolated groups of tsetse.

The study sites spanned a zone of contact between two divergent groups of G. fuscipes co-occuring in Uganda [9,10]. Sites MK and OT were situated north of the zone of contact and flies here were expected to possess solely northern mtDNA haplotypes. Sites BN, JN, and MS were located at the zone of contact and flies here were expected to possess both northern and southern haplotypes. Sites BU and OK were located south of the zone of contact and flies at these sites were expected to possess exclusively southern mtDNA haplotypes.

DNA Extraction

DNA was extracted from tsetse legs using NucleoSpin 96 Tissue Kits (Clontech, Mountain View, CA) or DNeasy kits (Qiagen, Valencia, CA) following the manufacturer's protocols.

Mitochondrial DNA sequencing

PCR was used to amplify a 570 bp fragment of mtDNA from a random subset of flies from each population using the primers COIF1 (CCT CAA CAC TTT TTA GGT TTA G) and COIIR1 (GGT TCT CTA ATT TCA TCA AGT A) as described by [10]. We amplified COIF1/COIIR1 in a 25 μl PCR reaction containing 1 × buffer (GoTaq colorless, Promega), 0.8 mM each dNTP, 0.4 mM primers, 1.5 mM MgCl2 and 0.5 U Go Taq polymerase. The amplification involved a denaturation step at 95°C for 8 min, followed by 50 cycles each at 94°C for 30 s, 51°C for 30 s, 72°C for 45 s, with a final extension step at 72°C for 7 min. PCR products were sequenced using an ABI Model 3730 automated sequencer (Applied Biosystems, Foster City, CA, USA). Electropherograms were visually inspected and sequences were trimmed to remove poor quality data. The resulting sequences (530 bp) were aligned by eye using the computer program Sequencher 4.2.2 (Gene Codes Corporation).

Microsatellite genotyping

We genotyped individual flies at 16 loci. We used 11 of the 13 loci described by [10], excluding D05 and Pgp17 due to possible null allele problems. We also employed five new dinucleotide loci identified in the G. morsitans genome and optimized for use in G. fuscipes: GmmA06, GmmB20, GmmD15, GmmL03, GmmL11 [14]. Amplifications were performed with fluorescently labeled forward primers (6-FAM, HEX and NED) using a touchdown PCR (10 cycles of annealing at progressively lower temperatures from 60°C to 51°C followed by 35 cycles at 50°C) in 12.5 μl reaction volumes employing 1 × buffer, 0.8 mM dNTPs, 2.0 mM MgCl2 and 0.5 U Go Taq polymerase. PCR products were multiplexed in groups of two or three and then genotyped on the ABI 3730 automated sequencer. Alleles were scored using the program Genemarker v 5.0 (SoftGenetics) with manual editing of the automatically scored peaks.

Marker validation and genetic diversity

Microsatellite loci were evaluated for Hardy Weinberg equilibrium (HWE) and linkage disequilibrium (LD) using Genepop version 4.0 [15]. Markov chain parameters were set at 10,000 dememorizations, 1000 batches, 10,000 iterations per batch for HWE and 100,000 dememorizations,1000 batches,10,000 iterations per batch for LD. We used the method of [16] as implemented in MultiTest v.1.2 to correct for multiple tests. Locus- and population-specific estimates of microsatellite allele frequencies were generated using the program Genalex version 6.3 [17]. We used the program FSTAT version 3.1 [18] to calculate allelic richness and the program Arlequin v. 3.5 [19] to calculate observed (Ho) and expected (He) heterozygosity for each population. DnaSP version 5.0 [20] was used to calculate mtDNA haplotype diversity (Hd) and nucleotide diversity (π).

Temporal genetic differentiation and population stability

For microsatellite data, we used Jost's DEST [21] to quantify genetic differentiation between populations and between temporal samples from the same population. DEST provides a less-biased estimate of differentiation than FST and related statistics, especially when estimated using highly polymorphic microsatellite loci [22]. Locus-specific calculations of DEST were performed using the web-based program SMOGD [23] and then averaged across loci. For mtDNA data, we used Fisher's exact test and the statistical software SAS version 9.1 to test for differences in haplotype frequencies among temporal samples from the same population. For both microsatellite and mtDNA data, we performed an analysis of molecular variance (AMOVA) as implemented in Arlequin v. 3.5 [19] to characterize the proportion of the variance in microsatellite allele frequencies or haplotype frequencies that was attributable to differences in date of sampling. For this analysis, we used only the two samples from each population that were separated by the longest time interval.

We estimated current effective population sizes based on temporal changes in microsatellite allele at all seven sites. The effective size of a population (Ne) is defined as the size of an ideal population (i.e., one of constant size, discrete generations, and negligible selection and gene flow), which would exhibit the same genetic characteristics as the population at hand [24]. Ne, therefore, reflects the rate of change in gene frequencies due to random drift alone [25]. We used two methods: the moment-based approach [26] and a likelihood approach implemented in TM3 [27]. Estimates were generated using the software NeEstimator [28]. For TM3, we employed 100,000 updates and a maximum Ne of 20,000.

For all analyses, we assumed that G. fuscipes undergoes approximately 8 generations per year using observations from colony flies (~7.3 generations per year, [29] ~8.5 generations per year at 25°C, [30] and those reported in other studies of the palpalis group (G. palpalis gambiensis and G. palpalis palpalis) in Guinea and Equatorial Guinea [5,31]. All populations were evaluated at an interval of at least one year (~8 generations apart). For four populations (BU, OK, BN, MK), we generated estimates at two different sampling intervals (0 to 8 generations, and 0 to 12 generations).

For each temporal sample in all seven populations, we also tested for an excess of heterozygosity relative to observed allelic diversity, which may be indicative of a recent bottleneck [32]. For each temporal sample, tests of heterozygosity excess were performed separately for each microsatellite locus. Significance was assessed across loci using Wilcoxon's test, which is the most appropriate test given the number of microsatellite loci evaluated. All tests were performed using the program BOTTLENECK [33].

Results

Marker validation and diversity

We genotyped a total of 667 tsetse flies at 16 microsatellite loci. We detected 17 values of FIS (out of 288) that exhibited significant departures from HWE at p < 0.05 (Additional file 1: Table S1). Assessed by locus, the number of significant FIS values observed was consistent with chance at an overall value of p < 0.05. Assessed by population, the number of significant FIS values observed was consistent with chance for all populations except the sample representing generation 12 from BN. Following sequential Bonferroni correction, only one locus pair exhibited significant linkage, and only in one population, confirming previous work showing that these loci were unlinked [10,14].

Microsatellite diversity was lowest in the samples from Mukongoro (MK) and highest in the samples from Bunghazi (BN). Allelic richness ranged from 3.0 to 4.4 and expected heterozygosity (HE) ranged from 0.418 to 0.609, (Table 1). MtDNA haplotype diversity was relatively low across samples with the number of haplotypes ranging from 1 to 4. As expected, nucleotide diversity was generally higher in populations from the zone of contact which were composed of flies with both northern and southern ancestry. We detected only two haplotypes that had not been previously reported [10]. Both of these haplotypes were recovered in population OT and differed by just one substitution from Hap26 or Hap27 [10].

Temporal variation in genetic diversity

Variation in allele frequencies by population and locus are depicted in Figure 2. Genetic differentiation between samples taken from the same population but at different times was extremely low, and uniformly lower than the differentiation observed between populations. DEST averaged 0.001 for temporally-spaced samples within populations, compared to 0.308 between populations (Additional file 2: Table S2). The largest values of DEST among temporally-spaced samples were observed in Masindi (MS generation 0 vs.13, DEST = 0.019 ± 0.022) and Otuboi (OT generation 0 vs. 11, DEST = 0.013 ± 0.007).

Figure 2

Microsatellite allele frequencies observed in seven populations of . Numbers after location code indicate the time interval (in generations) since the first sampling.

Mitochondrial haplotype frequencies also exhibited little change over time (Figure 3). We observed a significant change in haplotype frequencies only between the two temporally spaced samples from Junda (JN, p = 0.046). This was attributable to the loss of the two least common haplotypes in the sample representing generation 13.

Figure 3

Mitochondrial haplotype frequencies observed in seven populations of . Numbers after location codes indicate the time interval (in generations) since the first sampling. Only temporal samples from Junda (JN) differed significantly.

An analysis of molecular variance using both microsatellite allele frequencies and haplotype frequencies suggested that differences between temporally-spaced samples explained an insignificant amount of the overall genetic variation (Table 2). Differences among sites, on the other hand, contributed significantly to the overall variation in genetic diversity. The percent variation explained was greater in the case of mtDNA data, compared to microsatellite data.

Table 2

Results of an AMOVA testing for temporal genetic structure in seven populations of G. fuscipes sampled in 2008 and also in 2009/2010.

	d.f	Sum of squares	Variance components	% Variation	p
mtDNA 2008 vs. 2009/2010
Among temporal groups	1	0.4	-0.24319	-11.1	0.997
Among sites within groups	12	374.6	1.66812	76.4	0
Within sites	242	183.3	0.75752	35.7	0
microsatellites 2008 vs. 2009/2010
Among temporal groups	1	11.6	-0.15553	-3.2	1
Among sites within groups	12	986.2	1.07307	22.0	0
Within sites	1022	4057.4	3.96592	81.2	0

Results for microsatellites represent the weighted average over 16 loci. Variance components for which the expected value is zero may be slightly negative by chance.

Effective size

Estimates of Ne were generated for the seven populations based on microsatellite allele frequency changes observed among samples collected at different times from the same population. Momentum-based estimates ranged from 216 to infinity, but only the estimate from OT was bounded by a 95% confidence interval that did not include infinity (Table 3). Likelihood estimates ranged from 152 to 19,550 and all estimates were bounded by 95% confidence intervals that included 20,000, the maximum value of Ne considered (Table 3).

Table 3

Effective population size (Ne) and 95% confidence intervals (CI) for G. f. fuscipes populations.

Population	Interval sampled (generations)	Ne - moment	95% CI	Ne - likelihood	95% CI
BN	12	91926	575-infinity	13776	404-20,000
	8	1774	257.8-infinity	922.1	0-20,000
BU	12	7144	452-infinity	1203	0-20,000
	8	711	200-infinity	852	0-20,000
JN	13	infinity	293-infinity	2118	0-20,000
MK	12	1061	194-infinity	8024	393-20,000
	8	2170	161-infinity	11791	128-20,000
MS	13	443	154-infinity	19550	460-20,000
OK	12	infinity	807-infinity	19139	878-20,000
	8	infinity	378-infinity	18098	425-20,000
OT	11	216	117-439	152	95-20,000

Ne was calculated using the moments based temporal method of Waples (1989) and the likelihood approach of Berthier et al. (2002). For four populations, Ne was calculated using samples collected at intervals of both 8 and 12 generations from the initial sampling.

For populations MK and OK, estimates of Ne were similar regardless of whether the calculations were performed using data for generations 0 and 8 or generations 0 and 12. In populations BN and BU, however, estimates of Ne derived from the moment method differed by an order of magnitude depending on whether the sample representing 8 generations or 12 generations was included. In population BN, the estimate of Ne generated by the Likelihood method was similarly unstable.

Population bottlenecks

Following Bonferonni correction, seven samples (taken from Bunghazi (BN), Masindi (MS), Okame (OK) and Otuboi (OT)) exhibited significant signatures of a recent population bottleneck under the infinite allele model (IAM) model. Only one of these samples (MS generation 0) also tested positive for a bottleneck under the two phase model (TPM; Table 4). Samples from Busime (BU) and Mukongoro (MK) exhibited the least evidence for past bottlenecks (all but one p-value >> 0.05), however power to detect a bottleneck may have been low in MK on account of relatively low genetic diversity (Table 1).

Table 4

Significance of tests for population bottlenecks assessed using a Wilcoxon test under an infinite allele (IAM) or two-phase (TPM) model of microsatellite evolution.

Sample	p (IAM)	p (TPM)
BN - 0	0.00067	0.07193
BN - 8	0.00031	0.09686
BN - 12	0.00038	0.10388
BU - 0	0.22714	0.66061
BU - 8	0.17957	0.73776
BU - 12	0.07571	0.51102
JN - 0	0.02063	0.12619
JN - 13	0.06027	0.1514
MK - 0	0.0365	0.2106
MK - 8	0.16513	0.31934
MK - 12	0.2106	0.53296
MS - 0	0.00001	0.00258
MS - 13	0.01248	0.39098
OK - 0	0.01077	0.0535
OK - 8	0.00168	0.02094
OK - 12	0.00043	0.02899
OT - 0	0.00655	0.20187
OT - 11	0.00002	0.02396

Values in bold remained significant after Bonferroni correction (p < 0.0028).

Discussion

We assessed changes in genetic composition of seven tsetse populations in southeast Uganda in order to gain insight into the population dynamics of G. f. fuscipes. In general, our results provide evidence for temporal stability of G. f. fuscipes populations over the one to two year period that we examined. With the exception of just one or two populations discussed below, mitochondrial haplotype frequencies and microsatellite allele frequencies exhibited little change over time and effective population sizes were generally large.

Compared to other riverine species of tsetse, estimates of Ne for G. f. fuscipes were similar to or larger than estimates for G. palpalis palpalis in Equatorial Guinea [31] and 2 to 3 orders of magnitudes larger than estimates for G. p. gambiensis on islands off the coast of Guinea [5]. Values of Ne for G. f. fuscipes populations were also generally larger than estimates for a savannah species, G. pallidipes, in Kenya [34]. The large effective population sizes and overall stability of G. f. fuscipes populations support the hypothesis [35] that seasonal variation in tsetse numbers, in which larva develop in utero, should be relatively small, since they do not depend on surface water or moist media for breeding. Nonetheless, the lack of variation in genetic structure over time is surprising given the reduced abundance of G. f. fuscipes observed during the dry season in Uganda [12]. To reconcile our results with this observation, which may reflect the low efficiency of traps used for monitoring [36,37], we suggest that populations of G. f. fuscipes in dry season refugia remain large, and that seasonal invasion of marginal wet-season habitat (e.g., at Mukongoro, Bunghazi) must occur in waves of tsetse that are large enough to be representative of the refugia population. Large populations of pupa, which develop in the ground over a period of weeks, may also help to ensure the continuity of tsetse populations and would contribute to reducing the variance in genetic changes over time.

In contrast to the other populations, estimates of Ne were low for populations MS and especially OT, where both moment and likelihood methods produced values of only about 200. These values could be indicative of small populations. Ne may also be influenced by overlapping generations and temporal variance in reproductive success as well as the forces of selection, mutation and migration. In this study, however, the low values of Ne observed in these populations probably reflected small differences in the location of trapping sites used for the two temporal samples. Generation 13 from MS was sampled at a distance of about 4 km from the original site at which generation 0 was sampled. Likewise, generation 11 from OT was sampled at a single site that was 11-20 km from the relatively widely dispersed sites from which generation 0 was sampled. Thus, for these sites, which were the only two sites sampled at different locations across years, fine-scale spatial genetic variation could be responsible for the apparent temporal variation in gene frequencies, thus depressing estimates of Ne.

Given that genetic variation in MS and OT samples can probably be attributed to microgeographic variation, the change in genetic composition of the population at JN likely reflects the only significant temporal change observed in this study. Although microsatellite allele frequencies were largely invariate, mtDNA haplotype frequencies here differed significantly between generation 0 and generation 13. Junda (JN), along with sites BN and MS, lies along a narrow zone of contact between two long-diverged and historically-isolated groups of G. f. fuscipes [9,10]. In 2008, populations at all three of these sites harbored both "southern" and "northern" mtDNA haplotypes. Interestingly, in Junda, individuals with the "southern" haplotypes disappeared from the sample after 13 generations. This could be due to a particularly small population of females and stochastic variation in female reproductive success, although in tsetse, the latter is more likely to be true among males than females [5]. Mating success can also be influenced by Wolbachia, a symbiont that may impose mating barriers due to cytoplasmic incompatibility between infected and uninfected tsetse individuals [38], thereby biasing mating in favor of infected females and potentially producing mitochondrial sweeps [39]. Given the change in mtDNA observed at Junda, flies here should be examined for Wolbachia. If present, the zone of contact in Uganda may provide a unique opportunity to monitor symbiont-induced population changes over time.

Competing interests

The authors declare that they have no competing interests.

Authors' contributions

RE and JSB contributed to field collections, performed lab work, analyzed the data, and drafted an initial version of the manuscript. CH also performed lab work and data analysis. LMO, SA and AC helped design the study, coordinated fieldwork and provided guidance on the manuscript. All authors read and approved the final manuscript.

Additional file 1

Table S1. FIS values for the 16 microsatellite loci. Significance was assessed at p < 0.05 (*) and a Bonferroni-corrected value of p < 0.0028 (bold). Low variability precluded calculation of FIS in some populations (n/a).

Click here for file

Additional file 2

Table S2. Pairwise estimates of genetic differentiation (Jost's DEST) between samples taken from seven populations of G. f. fuscipes. Estimates of differentiation (below diagonal) and associated standard error (above diagonal) between populations of flies sampled at the same site but different times are shaded in grey.

Click here for file