Study area

Fieldwork was conducted in Estação Ecológica de Águas Emendadas (ESECAE), which is a 10547 ha Cerrado nature reserve close to Brasília, Brazil (15°32′–15°38′S and 47°33′–47°37′W). Mean annual precipitation in Brasília is 1·552 mm, with most of the rainfall occurring from November to January, while the dry season occurs during winter, which is from June to August (Maia and Baptista, Reference Maia, Baptista and Fonseca2008). Birds were captured in a 100 ha (1 km × 1 km) plot, which was delimited inside ESECAE. The plot area comprises 4 different vegetation types of the Cerrado biome: ‘cerrado’ (sensu stricto), ‘campo cerrado’ (open cerrado), ‘campo sujo’ (grassland with sparse shrubs) and ‘campo limpo’ (grassland without shrubs).

Field methods

Birds were captured between 2005 and 2009 in the reproductive season for most host species (August to December). Bird ringing was carried out periodically (25 sessions a year) by the members of the Laboratório de Ecologia e Conservação de Aves from Universidade de Brasília, as part of a long-term study on the demography and conservation of birds from the Cerrados of Central Brazil. Information on species reproductive and life-history parameters such as nest type (closed or opened), nest height, social system (pairs or groups), length of incubation, migratory status (resident or migratory species), and individual body mass were obtained from this long-term project database that consists of approximately 8500 ringed individuals and 3800 monitored nests between 2002 and 2010 (Marini, unpublished data).

Birds were captured in 10 to 20 mist-nets (mesh = 35 mm) of 10 or 12 m length and 2·5 m height. Nets were opened from sunrise until midday and were monitored every 30 min. Captured individuals were banded with metal rings issued by Centro de Pesquisas para Conservação de Aves Silvestres (CEMAVE/IBAMA), weighed with Pesola spring balance to the nearest 0·1 g and released. Species nomenclature follows the South American Classification Committee (Remsen et al. Reference Remsen, Cadena, Jaramillo, Nores, Pacheco, Pérez-Emán, Robbins, Stiles, Stotz and Zimmer2011). Approximately 50 μl of blood were obtained from the femoral or brachial vein using disposable needles or lancets and were stored in 1 ml of absolute ethanol (license issued by IBAMA, number 15256-1).

Molecular analysis

DNA samples were screened for Plasmodium and Haemoproteus infections using the primers 343F and 496R following Fallon et al. (Reference Fallon, Ricklefs, Swanson and Bermingham2003) with slight modifications. This protocol uses primers designed to amplify a 154 nucleotide segment of RNA-coding mitochondrial DNA (Fallon et al. Reference Fallon, Ricklefs, Swanson and Bermingham2003). Since Leucocytozoon, another closely related blood parasite, was absent from this bird community in a previous study (Fecchio et al. Reference Fecchio, Lima, Silveira, Braga and Marini2011), this genus of blood parasites was not screened for. PCR reactions were run with TaKaRa™ (TaKaRa Bio USA, Madison, WI, USA) reagents in 10 μl volumes that contained the following final concentrations: 1X buffer, 2 mmol MgCl₂, 0·25 mmol of each dNTP, 0·4 μmol of each primer, and 0·25 units Taq polymerase. One μl of the extracted DNA was used in each reaction. Thermal cycling conditions were as follows: initial denaturation of 2 min at 94 °C followed by 35 cycles of 50 sec denaturation at 94 °C, 50 sec annealing at 55 °C, and 25 sec extension at 72 °C. The reaction was completed with a 3 min extension at 72 °C. One μl of PCR product was run out on a 1% agarose gel stained with ethidium bromide and visualized with a UV light source.

To identify parasite infections to their cytochrome b haplotype, we first attempted to amplify a 702 bp segment of the cytochrome b gene for all samples determined to be infected by the 343F-496R screening procedure using the following primers and thermal cycling conditions: 3932F (5′-GGG TTA TGT ATT ACC TTG GGG TC- 3′; complement of 3932R in Olival et al. Reference Olival, Stiner and Perkins2007) and DW4 (Perkins and Schall, Reference Perkins and Schall2002), 1 μl of genomic DNA was subjected to an initial denaturation of 4 min at 94 °C, followed by 35 cycles of 94 °C for 20 sec, 49 °C for 10 sec, and 68 °C for 45 sec, and a final extension at 68 °C for 3 min. Reactions were run with TaKaRa™ reagents in 10 μl volumes that contained the final concentrations: 1X buffer, 0·2 mmol dNTP, 2 mmol MgCl₂, 0·04% BSA, 0·2 μmol of each primer, and 0·25 units of Taq.

A nested PCR was run using a 0·5 μl aliquot from the 3932F-DW4 reaction for those samples with weak or no bands, using the following primers and thermal cycling conditions: 413F and 926R (Ricklefs et al. Reference Ricklefs, Swanson, Fallon, Martinez-Abrain, Scheuerlein, Gray and Latta2005) under initial denaturation of 94 °C for 1 min and 28 cycles of 94 °C for 20 sec, 52 °C for 10 sec, and 68 °C for 50 sec, and a final extension at 68 °C for 7 min. Reactions were run with TaKaRa™ reagents in 20 μl volumes with the same concentrations as in the initial PCR reaction. Sequencing of the cytochrome b gene segments was carried out on an ABI 3100 Genetic Analyzer (Applied Biosystems, Foster City, CA, USA) according to the manufacturer's protocol. By using nested PCR our sequencing success was 90%; when lineage identity was unclear because of multiple infection, the sequence was removed from further analyses.

Phylogenetic analysis

Cytochrome b haplotypes (up to 663 bp long) sequenced in this study were aligned using Clustal X v2·0 (Larkin et al. Reference Larkin, Blackshields, Brown, Chenna, McGettigan, McWilliam, Valentin, Wallace, Wilm, Lopez, Thompson, Gibson and Higgins2007). We estimated a maximum likelihood phylogeny and assessed confidence in our clades by 100 bootstraps using the RAxML Black Box (Stamatakis et al. Reference Stamatakis, Hoover and Rougemont2008). RAxML Black Box uses as a default a GTR + G (gamma) model of evolution. Our tree was rooted with 6 mammalian Plasmodium species, that constitute the most appropriate outgroup to avian and lizard haemosporidians (Plasmodium spp., Haemoproteus spp., and Leucocytozoon spp.) based on outgroup-free rooting methods (Outlaw and Ricklefs, Reference Outlaw and Ricklefs2011). In addition to the RAxML analysis, we performed a Bayesian analysis in MrBayes v. 3.1.2 (Ronquist and Huelsenbeck, Reference Ronquist and Huelsenbeck2003) using the same evolutionary model (GTR + G) as for the ML analysis. We ran 3 million generations, sampling at every 500 steps, and discarded 50 000 generations as burn-in.

Defining haemosporidian lineages

Avian haemosporidian lineages were defined by comparing ca 600 nucleotides from the cytochrome b gene of the parasite. Because most of the lineages obtained in this study are not described species, they were defined, primarily, in relation to genetic variation and host species (Fallon et al. Reference Fallon, Bermingham and Ricklefs2005). There is no general agreement on the level of genetic variation of mitochondrial DNA required to define a species of avian haemosporidia (Ricklefs et al. Reference Ricklefs, Fallon and Bermingham2004, Reference Ricklefs, Swanson, Fallon, Martinez-Abrain, Scheuerlein, Gray and Latta2005). For instance, haemosporidia in birds with only 0·6% difference in their nucleotide sequences have been considered different lineages (Ricklefs et al. Reference Ricklefs, Swanson, Fallon, Martinez-Abrain, Scheuerlein, Gray and Latta2005), while up to 3·1% differences have been noted for cryptic Plasmodium species in Caribbean Anolis lizards (Perkins, Reference Perkins2000). In mammals, less than 1% divergence in cytochrome b sequences was observed in different Plasmodium species (Escalante et al. Reference Escalante, Freeland, Collins and Lal1998), and lineages differing by only 0·5% were shown to be evolutionarily independent in another bird study (Bensch et al. Reference Bensch, Pérez-Tris, Waldenström and Hellgren2004). In the present study, sequences differing in at least 3 nucleotides (0·5%) were considered distinct lineages.

Host specificity

To estimate host specificity we calculated the S ^*_TD index proposed by Poulin and Mouillot (Reference Poulin and Mouillot2005)

$$S_{TD}^* = \displaystyle{{\sum {\sum\limits_{i \lt j} {\omega _{ij} (\,p_i p_j )}}} \over {\sum {\sum\limits_{i \lt j} {(\,p_i p_j )}}}} $$

where the summations are over the set {i = 1, … , s; j = 1, … , s, such that i < j and s is the number of host species used by the parasite}, ω _ij is the pairwise genetic distance between host species i and j, p _i and p _j is the prevalence of the parasite lineage on host i and j, respectively. S ^*_TD is the weighted average of the genetic distance between hosts, regardless of the number of host species parasitized. In place of pairwise genetic distance we used Linnean taxonomic ranks, dividing the taxonomy into orders, suborders, families, subfamilies, genera, and species. For example, 2 species in the same genus have a distance of 1, whereas 2 species in the same order, but different suborders, have a distance of 5. The maximum value possible is 6, representing 2 species from different orders. The index ranges from a value of zero (host species specialist) to 6 (evenly distributed host order generalist). Lineages recovered only once receive a zero by default; however, with fewer than 2 samples it is impossible to assess host specialization. Moreover, high levels of uncertainty of the S ^*_TD value exists in poorly sampled lineages both because host range (estimated as the number of infected host species) is strongly related to sample size of haemosporidian lineages (Ricklefs et al. Reference Ricklefs, Swanson, Fallon, Martinez-Abrain, Scheuerlein, Gray and Latta2005) – it is unlikely that we have recovered poorly sampled lineages from all compatible hosts – and because we do not have an accurate estimate of prevalence in small samples. Additionally, some host species are poorly sampled, meaning that we are likely to miss rare haemosporidian lineages infecting these hosts (Ventim et al. Reference Ventim, Morais, Pardal, Mendes, Ramos and Pérez-Tris2012). However, the S ^*_TD gives a more biologically meaningful approximation of a lineage's host breadth than simply adding up the number of host species it infects.

Statistical analyses

To determine the influence of bird species traits on blood parasite prevalence, Generalized Linear Mixed-effect Models (GLMM) were fitted to the data using the ‘lmer’ function implemented in the R package lme4 (Bates and Maechler, Reference Bates and Maechler2010). For each GLMM (one for Plasmodium and one for Haemoproteus (Parahaemoproteus)), the response variable was the infection prevalence for each of the 14 most abundant bird species (minimum sample size = 14 individuals), which was inserted in R as a 2-vector response variable that included the number of infected individuals as well as the number of non-infected individuals (Crawley, Reference Crawley2007). We chose this cut-off based on Jovani and Tella (Reference Jovani and Tella2006). These authors argue that around 15 individuals per species could be used as a reasonable trade-off between not losing too much data from analyses and maintaining acceptable levels of uncertainty. Each model was assigned a binomial error distribution and a logit link function. To control for the lack of statistical independence in the data because of phylogenic effects, variance was decomposed across species by using a hierarchical taxonomic effect (Family/Genus) that was inserted in the models as a random term, while species biological parameters (nest height, species weight, social system, nest type, incubation time and migratory status) were considered as fixed effects. Because of the low number of species that were either closed or cavity nesters, these 2 types of nests were grouped into 1 single category.

Before models were fitted, the numeric variables were standardized by subtracting the mean and dividing by 2 S.D., except for the categorical variables, thus, effect sizes of models are on a comparable scale (Gelman, Reference Gelman2008). Data exploration was done in accordance to the protocol of Zurr et al. (Reference Zurr, Ieno and Elphick2010) to check for possible collinearity among the variables, homogeneity of variance, and outliers before models were fitted. GLMM assumptions were verified in the global model (because we used AIC_c for model selection, see below), which included all the variables without interaction, by checking diagnostic plots of the residual distribution and conditional modes for the random effects. Estimated model variance was verified for overdispersion; including the random taxonomic effect reduced overdispersion (i.e. when compared with a GLM), and neither of the GLMM were overdispersed. Our models were consistent with mixed model assumptions.

Several explanations have been proposed for the observed variation among host species in parasite prevalence (see Introduction section) and we have considered many of these hypotheses in the set of candidate models to be tested, including a null model (no fixed effects, only the intercept with random effects) and a global model (all fixed factors in addition to the random effects). We did not include any models with interactions among the 6 fixed factors because we do not know of any empirical support or theoretical basis a priori to include their interactions. However, we did include all possible variable associations without interaction terms in the set of candidate models (62 models) that together with the null and global model gave a total of 64 models.

We selected models using second-order Akaike's Information Criteria (AIC_c), which is derived from information theory and measures the relative strength of evidence for each alternative model in a set of candidate models (Burnham and Anderson, Reference Burnham and Anderson2002). AIC_c measures the relative goodness of fit of a model, while taking the number of parameters in each of the models into consideration, and models are penalized with the addition of parameters. Thus, AIC_c selects a model that has a good fit with a minimum number of parameters (Mazerolle, Reference Mazerolle2004).

Model selection and model averaging were done using the R package AICcmodavg (Mazerolle, Reference Mazerolle2010) and all other statistical tests were done in R (R Development Core Team, 2008).

Models were ranked using ΔAIC_c and Akaike's weight (w _i). Models with ΔAIC_c values up to 2 are generally considered to present substantial support, but models with ΔAIC_c up to 7 may contain sufficient information to be considered plausible hypotheses (Burnham et al. Reference Burnham, Anderson and Huyvaert2011), while w _i represents the probability of a model being the best among the candidate models and varies from 0 (no support) to 1 (complete support) (Burnham and Anderson, Reference Burnham and Anderson2002; Mazerolle, Reference Mazerolle2004). Multi-model inference was carried out in the a priori set of models that presented ΔAIC_c equal or lower than 7 (Burnham et al. Reference Burnham, Anderson and Huyvaert2011). Estimates of parameters by model averaging and calculations of the unconditional standard errors and unconditional confidence intervals (C.I.) were done as described by Burnham and Anderson (Reference Burnham and Anderson2002).