Philornis downsi (family Muscidae; subfamily Azeliinae;
tribe Reinwardtiini) is a semi-haematophagous obligate avian parasite
in its three larval stages, whereas adult flies are non-parasitic and
feed on organic matter .
Adults lay eggs inside the nares of newly hatched nestlings (usually at
one to three days old), which hatch into first instar larvae [17,24]. Second and third instar larvae attach externally and feed on nestling blood and tissues over four to six days . Most larvae of P. downsi appear
to reach their third instar phase at the time of host fledging. The
larvae pupariate at the base of the nesting material and remain for
approximately two weeks before emerging as adult flies [13,25].
Study area and sample collection
Philornis downsi were collected from three islands of the Galápagos: Santa Cruz (986 km2; 0° 37'S, 90° 21'W), Floreana (173 km2, 1° 28'S, 90° 48'W), and Isabela (4588 km2,
0° 58'S, 90° 58'W). Fly samples were collected from nests during the
January to March finch breeding season in 2004, 2005 and 2006 from two
contrasting habitats, the arid lowlands (0–100 m asl) and the humid
highlands (300–600 m asl) (Table 1) [see also [26,27]]. The lowlands are characterised by low rainfall, and are dominated by the trees Acacia macracantha, Bursera graveolens, Croton scouleri, Opuntia spp., Pisonia floribunda, and Zanthoxylum fagara . In contrast, the highlands have much higher rainfall [28,29], abundant moss and lichen, and are dominated by the endemic tree Scalesia pedunculata, or S. cordata (Asteraceae) on Isabela Island.
Table 1. Sample sizes of bird nests and P. downsi individuals.
We sampled from one site in each habitat on both
Floreana (lowlands, adjacent to the town of Puerto Velasco Ibarra: 1°
16'S, 90° 29'W; highlands, base of Cerro Pajas: 1° 17'S, 090° 27'W)
and Isabela (lowlands: adjacent to town of Puerto Villamil: 0° 57'S,
91° 00'W; highlands: 0° 50'S, 91° 01'W), while on Santa Cruz we sampled
from three sites in the lowlands: (1) Garrapatero: 0° 39'S, 90° 28'W;
(2) Itabaca: 0° 29'S, 90° 17'W; (3) Punta Estrada, near Puerto Ayora: 0
° 45'S, 90° 18'W, and one site in the highlands (Los Gemelos: 0° 37'S,
90° 22'W) (Figure 1). All sample sites were approximately 2000–4000 m2, except for the highland site on Isabela, where our sample site was only 100 m2 because habitat fragmentation has reduced the Scalesia forest
to small remnant patches. The distance between highland and lowland
sites was much shorter on Floreana (3–5 km) than on Santa Cruz and
Isabela (both 15–25 km), while on Santa Cruz, the distance between all
four sites (1 highland, 3 lowland) varied between 15 and 27 km. Data
were obtained from all three islands in 2004, from just Floreana in
2005, and from Santa Cruz and Floreana in 2006 (Table 1).
Figure 1. Map of Santa Cruz, Floreana, and Isabela Islands with sampling locations.
Sampling sites on Santa Cruz: S1 = highlands; S2 = lowlands, Punta
Estrada; S3 = lowlands, Garrapatero; S4 = lowlands, Itabaca. On
Floreana and Isabela, one site each in the lowlands (L) and highlands
(H) are indicated.
For the purpose of our study, larvae, puparia and puparia cases were sampled from 64 bird nests of five Darwin finch species (Geospiza fuliginosa, n = 25, Geospiza fortis, n = 15, Camarhynchus parvulus, n = 3, Camarhynchus pauper, n = 4; Cactospiza pallida, n = 1), while one nest was opportunistically sampled from each of the Galápagos mockingbird (Nesomimus parvulus) and the yellow warbler (Dendroica petechia aureola). Fourteen recently fledged nests were sampled for P. downsi where
the finch species was unknown. GPS coordinates were recorded at each
nest location. Inactive nests were collected and sealed in individual
plastic bags and later dismantled for counting of P. downsi individuals. All flies were immediately preserved in 95% ethanol.
DNA extraction and microsatellite typing
DNA extraction was carried out using the salting out procedure described in  with the exception that all samples (3 mm2 tissue
from each individual) were homogenised and washed three times in 10 mm
TRIS prior to digestion with Proteinase K to remove traces of ethanol,
excess lipids, and other potential contaminants. Across all three
islands, 1012 P. downsi individuals (larvae and pupae) were genotyped (Table 1) using eight microsatellite markers : Pd1 [GenBank: EF608562] Pd2 [EF608556], Pd4 [EF608557], Pd6 [EF608564], Pd7 [EF608558], Pd8 [EF608555], Pd9 [EF608561], Pd10 [EF608563]. Multiplex PCR conditions were followed as described in Dudaniec et al. .
Samples were genotyped on an ABI 3730 capillary electrophoresis DNA
analyser (Applied Biosystems). A fluorescently labeled size standard
(GS500 (-250) LIZ) was run with the samples and alleles were scored
using GENEMAPPER version 3.7 (Applied Biosystems). To
minimise and estimate genotyping error, each run of the DNA analyser
contained eight repeated samples and a control sample run each time. In
total, this resulted in 70 individual samples (14.5% of all samples
genotyped) being re-amplified and genotyped at least once.
Mitochondrial DNA sequencing
An 822-bp region of the 3' end of the CO1 gene was amplified in five P. downsi individuals
collected from Santa Cruz (1 highlands), Floreana (1 highlands, 1
lowlands), and Isabela (1 highlands, 1 lowlands). Samples were
amplified using primers M202 (forwards, C1-J-1751 ) and M70 (reverse, UEA10 ). Amplifications were performed in 10× TaqGold buffer, 25 mM MgCl2, 10 mM total dNTP's, 200 nM each primer, 0.2 U TaqGold
polymerase, and 10–50 ng DNA. Amplification conditions were an initial
denaturation at 94°C for 9 min, followed by 34 cycles of 94°C for 45 s,
55°C for 45 s, 72°C for 1 min, with a final extension of 72°C for 6
min. Sequencing was performed using the ABI Prism™ Big Dye Terminator
Cycle sequencing kit (Applied Biosystems) according to the
manufacturer's instructions. Products were sequenced on ABI 3700
(version 3.7) automated DNA sequencers. SeqEd (version 1.0.3) (Applied
Biosystems) was used to edit chromatogram files to determine
bi-directional consensus sequences and to manually align sequences
Allele frequencies and data set construction
We calculated allele frequencies using RELATEDNESS 5.0.8 
by randomly selecting one individual per sample (n = 64) to eliminate
the possibility of including related individuals (a sample is defined
as all P. downsi individuals collected from a single bird's
nest). Exact tests were performed for each microsatellite locus to test
deviation from Hardy-Weinberg equilibrium using GENEPOP. All loci were in Hardy Weinberg equilibrium after sequential Bonferroni correction  and these allele frequencies were used for all further analyses. Genetic relatedness among P. downsi offspring
within nests of Darwin's finches is low, and the individuals found
within each nest are produced by up to approximately five ovipositing
females that have each mated with between one and five males (as found
by sib-ship reconstruction analysis by Dudaniec et al. in review). To
eliminate the effect of sibs in the data, we selected unrelated
individuals that were identified using the sib-ship reconstruction
method implemented in the program COLONY 1.2 . Each sample of P. downsi individuals taken from an independent bird nest was run in COLONY
1.2, which uses a maximum likelihood method that partitions individuals
into pure full-sib families (i.e. monogamous female parent), or
full-sib families nested within half-sib families (i.e. polyandrous
female parent) using progeny genotypes without known parental genotypes
Three runs were performed per sample with different random seed numbers
(12, 80, and 243) to ensure data convergence, and a conservative error
rate of 5% was implemented based on evidence from the re-genotyping of
70 individuals, in which genotyping error ranged from 0–5% across loci.
We selected one individual per reconstructed maternal family (i.e.
one family = the offspring assigned to one putative female parent). In
nested-half sib families (i.e. one mother, multiple fathers),
individuals were only selected from full sib families with the largest
number of members that had the highest posterior probability. Only
individuals genotyped at all eight loci were included in the analysis
and individuals were not sampled from families that contained Class I
or Class II typing errors (identified by COLONY 1.2) .
These criteria resulted in a sample size of 158 individuals sampled
from 63 bird nests (with between one and six unrelated individuals per
sample) (Table 1).
To examine the probability that two randomly selected individuals
from the same population will have the same multi-locus genotype, a
Probability of Identity (PI) analysis was performed using GIMLET.
The output is a cumulative multi-locus PI value, estimated both with
and without sample size correction. PI values were calculated for the
dataset of 158 individuals using equations of unbiased PI, which
assumes that individuals are unrelated, and PI for sibs, which assumes
that all individuals are siblings .
Inter-island genetic differentiation
Heterozygosity, and pairwise Fst  was calculated to examine genetic differentiation between islands (Santa Cruz, Floreana, Isabela) using MICROSATELLITE ANALYSER (MSA) 4.05 . Genotypic differentiation was tested between islands using option 3 with 10 0000 Markov chain iterations in GENEPOP. P-values for multiple tests were adjusted using sequential Bonferroni correction . The AMOVA method  was conducted in GENALEX version 6 
to partition the total genetic variation into three levels: among
islands, among individuals, and within individuals using the
Codom-genotypic distance calculation and 9999 permutations.
Population bottleneck analysis
Recently colonised species may experience a population bottleneck,
resulting in a reduction in the number of alleles and expected
heterozygosity at polymorphic loci. However, alleles may be lost at a
faster rate than the loss of heterozygosity, so observed heterozygosity
is higher than the expected heterozygosity at equilibrium . The program BOTTLENECK version 1.2.02  was used to test for the presence of a recent population bottleneck for P. downsi by
analysing within-population heterozygosity and allele frequency using
the constructed dataset of individuals sampled from all islands. Both
the stepwise mutation model (SMM) and two-phase model of mutation (TPM)
were used, with the latter model being considered the most appropriate
for microsatellites. The variance for the TPM was set at 5% and the
proportion of SMM in TPM was set at 95% .
To determine differences in gene diversity across loci, the Wilcoxin
sign-rank test was used as recommended for data sets with less than 20
loci (with 10 000 permutations) .
We also examined the allele frequency distribution in order to see
whether it is approximately L-shaped (as expected under mutation-drift
equilibrium) or not (indicating that a recent bottleneck has provoked a
mode shift) as described in .
Genetic structure among islands
For inferring genetic structure among the three sampled islands, we
conducted two complementary individual-based Bayesian clustering
analyses using STRUCTURE 2.1  and the landscape genetics program GENELAND without a priori knowledge
of population units and limits. Both software packages implement a
Bayesian clustering method that uses a MCMC technique to define the
number of populations in a sample that are at Hardy Weinberg
Equilibrium. The methods implemented in these two programs differ in
that GENELAND determines the optimal number of
populations or 'clusters' and then allocates individuals
(probabilistically) to these clusters using geographic coordinates,
whereas STRUCTURE carries out the allocation
sequentially for different numbers of clusters, and then flags the
number of clusters with the highest likelihood . In STRUCTURE,
the following run parameters were used: admixture without population
information used, correlated allele frequency model, a burn-in period
of 100 000 simulations followed by a run length of 1 million Markov
Chain Monte Carlo (MCMC) simulations and three iterations for each
number of potential clusters (defined as k = 1–5) to check for consistency of results. Estimation of k was taken to be the values of k with the highest Pr (X|k).
In contrast to STRUCTURE, the algorithm implemented in GENELAND is considered to be a powerful clustering method under conditions of low genetic differentiation among populations [51,52].
The model infers genetic discontinuities between populations in space
from multilocus genotypes obtained from geo-referenced individuals [49,53].
All individuals from the same sample (i.e. same bird nest) were
allocated the same GPS coordinates. GPS coordinates were available for
57/63 nests (Santa Cruz: n = 18; Floreana: n = 36; Isabela: n = 3) (138
individuals in total). Samples for which GPS coordinates were missing
were excluded from the analysis. To firstly infer the number of genetic
clusters (k) in our data set, we used the Dirichlet model,
which assumes independent allele frequencies with the following
parameters: 1000 000 MCMC iterations, uncertainty attached to spatial
coordinates = 0, variable number of populations = TRUE, minimum k = 1, maximum k = 5, and spatial information included in the model = TRUE. This procedure was performed three times to establish consistency of k across runs. The established k was
then run five times to check the consistency of individual assignment
to the inferred populations across runs. The same parameters were used
but k was fixed at the modal number found in the first
analyses. These five runs were post-processed (with a burn-in of 1000 ×
100 iterations) to obtain posterior probabilities of population
membership for each individual. Consistency of results across the five
runs was checked visually.
Inferred populations were further examined for heterozygosity, allelic richness (corrected for sample size), observed (Ho) and expected (He) heterozygosity, inbreeding coefficients (Fis), and genetic differentiation (estimated using Fst) using FSTAT v. 2.9.3 .