Seasonal dynamics of haemosporidian (Apicomplexa, Haemosporida) parasites in house sparrows Passer domesticus at four European sites: comparison between lineages and the importance of screening methods.

Infectious diseases often vary seasonally in a predictable manner, and seasonality may be responsible for geographical differences in prevalence. In temperate regions, vector-borne parasites such as malaria are expected to evolve lower virulence and a time-varying strategy to invest more in transmission when vectors are available. A previous model of seasonal variation of avian malaria described a double peak in prevalence of Plasmodium parasites in multiple hosts resulting from spring relapses and transmission to susceptible individuals in summer. However, this model was rejected by a study describing different patterns of seasonal variation of two Plasmodium spp. at the same site, with the double peak only apparent when these species were combined. Here, we assessed the seasonal variation in prevalence of haemosporidian parasites (Plasmodium, Haemoproteus and Leucocytozoon) in house sparrows (Passer domesticus) sampled across 1 year at four temperate European sites spanning a latitudinal range of 17°. We showed that parasite prevalence and diversity decreased with increasing latitude, but the parasite communities differed between sites, with only one Plasmodium lineage (P_SGS1) occurring at all sites. Moreover, the nested PCR method commonly used to detect and identify haemosporidian parasites strongly underestimated co-infections of Haemoproteus and Plasmodium, significantly biasing the pattern of seasonal variation, so additional molecular methods were used. Finally, we showed that: (i) seasonal variation in prevalence of haemosporidian parasites varied between study sites and parasite lineages/species/genera, describing further cases where the double peak model is not met; (ii) the seasonal dynamics of single lineages (P_SGS1) varied between sites; and (iii) unexpectedly, seasonality was greatest at the southernmost site, a pattern that was mostly driven by lineage H_PADOM05. Limitations of the genotyping methods and consequences of pooling (parasite lineages, sites and years) in studies of haemosporidian parasites are discussed and recommendations proposed, since these actions may obscure the patterns of prevalence and limit ecological inferences.

Introduction

Many infectious diseases are notoriously influenced by season, often showing cyclic, predictable patterns of incidence (Dowell, 2001, Altizer et al., 2006). These include human and wildlife pathogens with direct or airborne transmission such as poliomyelitis and influenza, as well as vector-borne diseases such as dengue, encephalitis and malaria (Randolf et al., 2000, Dowell, 2001, Hoshen and Morse, 2004, Li et al., 2019). Seasonal variation in disease dynamics has been associated with multiple drivers (reviewed by Altizer et al., 2006). Periodic variation in the density of hosts and in social interactions affect the transmission of many diseases, including avian influenza on the migration gatherings of birds (Krauss et al., 2010) and mycoplasmal conjunctivitis associated with flock size in house finches (Carpodacus mexicanus; Altizer et al., 2004). Disease resistance of hosts may also vary seasonally, for instance due to physiological changes (such as sex hormone production in adults during the breeding season) or due to differences in age compositions of the population caused by the seasonal production of immunologically naïve juveniles (Beaudoin et al., 1971, Dowell, 2001, Altizer et al., 2004). In addition, the predictable seasonal variation in weather is known to affect, for instance, availability, number and infectivity of vectors such as mosquitos and ticks (Randolf et al., 2000, Li et al., 2019). Importantly, seasonal drivers of disease dynamics can be responsible for geographical differences in transmission, timing and persistence of diseases (Altizer et al., 2006, Li et al., 2019), with more pronounced epidemics generally being positively associated with latitude (e.g. Cook et al., 1990, Dowell, 2001) or with specific measures of seasonality (e.g. annual amplitude of primary production; Lisovski et al., 2017; for an exception see Altizer et al., 2004).

Theory suggests that vector-borne pathogens in seasonal environments (where vectors are available for only a few months) should evolve lower virulence and transmission (Cornet et al., 2014), as well as a time-varying strategy so that greater investment in transmission occurs when vectors are available (Pigeault et al., 2018). Haemosporidian parasites of the genera Plasmodium, Haemoproteus and Leucocytozoon have persistant stages in their life cycle that can remain latent in various host tissues (Valkiūnas, 2005), thereby being largely protected from the host immune system during the winter months when vectors are generally not available and transmission cannot occur. In line with the theory, it was shown for the human pathogen Plasmodium vivax that the relapse periodicity was strongly dependent on latitude, being frequent in the tropics and prolonged in temperate regions (interpreted as an adaptation to seasonal changes in vector survival in order to optimise transmission; Battle et al., 2014).

Unlike human malaria, for which the timing and level of incidence can be reliably predicted using deterministic models incorporating weather factors (e.g. Hoshen and Morse 2004), relatively little is known about the epizoology of the huge diversity of avian haemosporidian parasites (for an exception see Samuel et al., 2011). Multiple studies have shown an increase in prevalence of Plasmodium parasites in the blood of resident birds during spring and summer in temperate regions (e.g. Micks, 1949; Janovy, 1966; Beaudoin et al., 1971, Cosgrove et al., 2008), although prevalence can also remain high for some parasites during winter (e.g. Dunn et al., 2014). Beaudoin et al. (1971) developed a functional model of seasonal variation of avian malaria for temperate regions that explained the occurrence of a double peak of prevalence of Plasmodium parasites in multiple bird hosts. The first peak of occurrence was thought to be caused by the relapse of infections that had been latent during the winter in adult birds (in tissues other than blood; Applegate, 1971, Beaudoin et al., 1971; see also Valkiūnas, 2005, Svensson-Coelho et al., 2016), and this peak could occur in early spring before most vectors become available (Janovy, 1966). Relapses were experimentally shown to be dependent on the physiology of the birds and particularly on their hormone levels (Applegate, 1970, Valkiūnas et al., 2004), although relapses could potentially also be stimulated by mosquito bites or be driven by an internal clock of the parasites (Pigeault et al., 2018). In addition, spring relapses were shown to increase the rate of spread of malaria to the mosquito vectors, thereby being able to transmit the disease to uninfected, susceptible birds (Applegate et al., 1971; see also Cornet et al., 2014). The second peak of prevalence, occurring during summer/autumn, was mostly dependent on the transmission of malaria parasites to immunologically naïve juveniles, which by then constituted a large proportion of the host populations (Beaudoin et al., 1971). The occurrence of this second peak was, however, much more variable between years and appeared to be affected by weather conditions and the consequent availability of vectors (although the parasites are able to persist as latent infections in the avian reservoir; Beaudoin et al., 1971). Various questions regarding this model have to a large degree been left unsolved, particularly the role of migratory bird species (which had higher Plasmodium prevalences) in introducing parasites into the bird community (Beaudoin et al., 1971).

Cosgrove et al. (2008) examined the seasonal variation of Plasmodium parasites in blue tits (Cyanistes caeruleus), showing a double peak of prevalence when the two commonest Plasmodium spp. (Plasmodium relictum and Plasmodium circumflexum) were combined. However, P. relictum had similar levels of prevalence throughout the year, whereas P. circumflexum showed double peaks (similar to those outlined by Beaudoin et al., 1971) only when both age classes of blue tits were combined. Therefore, although no fundamental functional aspect of Beaudoin’s model was challenged by Cosgrove et al. (2008), it was evident that the actual seasonal pattern of malaria prevalence in birds is variable between parasites, even at a single location within the temperate region and within a single host species. Additional studies have analysed seasonal variation of malaria parasites in birds (e.g. Atkinson and Samuel, 2010, Hellgren et al., 2013, Podmokła et al., 2014, Ishtiaq et al., 2017). However, several studies included host individuals belonging to different breeding populations (e.g. Hellgren et al., 2013; Pulgarín-R et al., 2019), thereby confounding seasonal and local effects (such as temperature and proximity to water; Wood et al., 2007; Lachish et al., 2011; Loiseau et al., 2013; Ferraguti et al., 2018). Moreover, combining samples from multiple years is a common practice, although this may potentially confound seasonal and annual variation (see Bensch et al., 2007).

In this study, we set out to investigate the seasonal variation of prevalence of haemosporidian parasites during a single year in a resident species, the house sparrow (Passer domesticus), studied in parallel at four locations spanning a latitudinal gradient of c.17° in the European temperate region. This is part of a larger study on local adaptation relative to parasite communities, and our aim was threefold. First, to test whether parasite prevalence and diversity decreased with increasing latitude, as has been shown previously in house sparrows, but using samples collected in multiple seasons and years (Marzal et al., 2011). Second, to evaluate whether seasonality in prevalence was in line with the seasonality at each location, assuming an increase with latitude at this geographical scale (see Lisovski et al., 2017). Third, to determine whether there were differences in seasonal patterns of occurrence between the commonest parasites of our study system. In addition, we explored and described limitations (that we came across) of the commonest molecular methods used to detect and identify haemosporidian parasites in birds, showing that these limitations can, in some circumstances, lead to false ecological inferences.

Materials and methods

Sample collection

Blood samples were collected from house sparrows in rural areas of Badajoz, Spain (38°83′20.64″N, −6°91′67.33″E), Burgas and Kalimok Biological Station, Bulgaria (43°85′67.37″N, 26°3′20.77″E and 44°01′15.82″N, 26°43′93.25″E, respectively), Trzcinica, Poland (49°7′25.89″N, 21°42′23.44″E and 51°24′29.74″N, 16°93′29.92″E), and at three sites in Skåne, Sweden (55°73′73.56″N, 13°6′29.87″E; 55°81′04.27″N, 13°58′81.21″E and 55°67′25.27″N, 13°53′82.46″E) located within a radius of 20 km. Birds were captured with mist nets and sampled every 2 months in each country between September 2016 and September 2017 inclusive; except in Bulgaria where no sampling was performed in September 2017. Additional samples collected at the same sites in June 2016 in Bulgaria (15 samples) and mid-August 2017 in Spain (11 samples) were used to compare the diversity of parasites and the molecular methods, but not in the assessment of temporal variation, as no comparable samples were available from other countries at the same time of the year.

Birds were individually marked with aluminium rings issued by the ringing centres of each country, sexed, aged (first year or adult), weighed and measured (wing length) whenever possible. Blood samples were also obtained from recaptured birds, ringed in previous sampling periods: 12 in Spain (out of 174 individuals), 27 in Bulgaria (n = 130), 17 in Poland (n = 150) and 22 in Sweden (n = 174). Blood samples were preserved in SET-buffer (0.15 M NaCl, 50 mM Tris, 1 mM EDTA) and kept at −40 °C until DNA extraction.

Genetic analyses

DNA was extracted from 250 µl of each preserved blood sample by using a standard ammonium acetate method; DNA pellets were dissolved in 100 µl of 1× TE Buffer (10 mM Tris, 1 mM EDTA) and its concentration was evaluated using a nano-drop spectrophotometer. Samples were then diluted to 10 ng/µl and stored at −40 °C.

The nested PCR described by Hellgren et al. (2004) was used to detect and identify blood parasites belonging to the genera Plasmodium, Haemoproteus and Leucocytozoon. The first PCR was performed in 25 µl reactions including 20–30 ng of template DNA, 1× PCR Buffer, 1.5 mM of MgCl2, 0.4 µM of each primer, 0.125 mM of each nucleotide, 0.5 units of Taq DNA-polymerase; and a thermal profile of 3 min at 94 °C, 20 cycles of 30 s at 94 °C, 30 s at 50 °C and 45 s at 72 °C, followed by 10 min at 72 °C. The second PCR included 1 µl of the product of the first PCR as well as the same quantities of the other reagents and only varied in that 35 cycles were performed instead of 20. Parasites were detected by running 2.5 µl the final PCR product on a 2% agarose gel for 25 min at 80 V, with infected individuals showing bands c.500 bp long. All 96 plates included negative controls (where water instead of DNA dilution was used for the first PCR), which never showed any bands in the gel.

PCR products in which parasites were detected were precipitated using 11 µl of ammonium acetate (NH4Ac, 8 M) and 37.5 µl of 95% ethanol, and dissolved in 5–25 µl of water (depending on the strength of the band visible on the gel). Precipitated PCR products were then sequenced with BigDye™ terminator cycle in 10 µl reactions including 0.75× Sequencing Buffer, 0.5 µM of the forward primer (HaemF), 1 µl of Big Dye® Terminator and 2 µl of the precipitated PCR product; and a thermal profile of 25 cycles of 10 s at 96 °C, 5 s at 50 °C and 4 min at 60 °C. The products of the sequencing reaction were precipitated by adding 2.5 µl of EDTA (125 mM) and 35 µl of 95% ethanol, and then sequenced on an ABI 3500 Genetic Sanger Sequencer (Applied Biosystems).

In addition, we used the multiplex PCR method described by Ciloglu et al. (2019) to assess whether co-infections of Haemoproteus and Plasmodium were being missed by the nested PCR method, which could have affected our results from Spain, where an Haemoproteus lineage proved to be very common. The multiplex method uses a combination of genera-specific primers in a single PCR, allowing us to detect whether Plasmodium, Haemoproteus and/or Leucocytozoon are present in the bird’s blood samples as they produce different-sized bands in the agarose gels. However, as the multiplex PCR method does not allow us to identify the specific lineages of parasites, new primers were designed to try to identify Plasmodium lineages that were present in co-infections with Haemoproteus (for which the nested PCR only identified Haemoproteus lineages). The Plasmodium-specific primers, which amplify the majority of the Plasmodium parasites infecting house sparrows in Europe, were designed using a whole cytochrome-b (cyt-b) sequence assembly of mitochondrial genomes of Plasmodium and Haemoproteus genera. Several primers located in dissimilar regions between the two genera but consistent among Plasmodium lineages were designed and tested. Finally, the primer combination Plas_4F (AATACCCTTCTATCCAAATCTTTTAA) and Plas_1R (TAAATAAACGACCATATAAAATGTAAATATC), which amplify a 368 bp fragment of cyt-b from Plasmodium parasites, was used in samples showing evidence of having co-infections of Plasmodium and Haemoproteus according to the multiplex PCR method. To do so, we performed PCRs in 25 µl reactions including: 50 ng of template DNA, 1× PCR Buffer, 1.5 mM of MgCl2, 0.4 μM of each primer, 0.125 mM of each nucleotide (dNTPs: deoxynucleotide triphosphate), 0.5 units of Taq DNA-polymerase. The thermal profile started with a short denaturation at 94 °C for 3 min, followed by 35 cycles at 94 °C for 30 s, 52 °C for 30 s and 72 °C for 45 s and finished with 72 °C for 10 min. PCR results were evaluated using 2.5 μl of each PCR product on a 2% agarose gel for 25 min at 80 V. All positive samples infected by Plasmodium parasites showed fluorescent bands (approx. 400 bp) and were subsequently precipitated and sequenced as described above.

DNA sequences were edited and aligned in Geneious© 11.1.5 (https://www.geneious.com), and compared with those available at the MalAvi database (Bensch et al., 2009). Mixed and co-infections (of the same or of different genera, respectively) were detected by the presence of double peaks in the electropherograms, and each parasite was identified by comparison with known parasite sequences described in MalAvi. A simple phylogenetic tree of all parasites detected in this study was produced in Geneious using the unweighted pair group method (UPGMA) procedure and the substitution model Tamura-Nei (which was the best of those available as determined by jModelTest 2.1.10; Guindon and Gascuel, 2003, Darriba et al., 2012). This tree allowed us to identify clades of lineages that have not been described morphologically but differed only by 1–2 bp from lineages belonging to known species; and such clades were considered species in our analyses (see Fig. 1).

Fig. 1

Unweighted pair group method tree of all the parasite lineages that were found in house sparrows, with bootstrap support (based on 100 replicates) of the major branches (performed in Geneious© 11.1.5). The lineages belonging to known species’ clades are indicated (MalAvi database, Bensch et al., 2009).

The Plasmodium-specific primers produced good sequences from all tested samples that were previously shown (using nested and multiplex PCRs) to have Plasmodium lineages (except for the lineage P_LINOLI01, which occurred in only one individual), whereas all tested samples infected only by Haemoproteus lineages were negative. However, due to limited overlap between the new sequenced fragment and the MalAvi fragment, this method only allowed us to distinguish the clades of Plasmodium spp., i.e. P. circumflexum, P. relictum, Plasmodium elongatum and Plasmodium cathemerium (see Supplementary Table S1), but not the specific MalAvi lineages. Whenever we failed to obtain sequences from putatively infected individuals, the sequences were too short or messy or the PCR methods produced inconsistent results, the whole procedure was repeated. New parasite lineages, and a few individuals for which the identification of the parasite remained doubtful, were also sequenced using the reverse primer (HaemR2). However, there remained some samples for which we could not identify the parasites. These samples (n = 3) were included in analyses of overall prevalence, but excluded from lineage-specific analyses.

Statistical analyses

All statistical analyses were performed in R v.3.5.1 (R Core Team, 2013. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL http://www.R-project.org/.). We used Generalized Linear Models (GLMs; with logit-link function and binomial distribution) in order to investigate whether the probability of infection (the binary dependent variable was defined as infected versus non-infected individual) differed between countries (Bulgaria, Poland, Spain and Sweden), sampling period (seven categories, from September 2016 to September 2017) and their interaction. Mixed models with ring number included as a random effect (to control for the non-independence of recaptured individuals) did not fully converge, but produced qualitatively similar results to GLMs including only captures, and so only first captures were analysed. As the probability of infection was not influenced by differences between sampling sites within Sweden (GLM Chisq test: P = 0.72) nor within Bulgaria (GLM Chisq test: P = 0.68), sites within countries were combined in all analyses.

We also used generalized additive models (GAMs), with binomial distribution and logit-link function, implemented in the R package mgcv (Wood, 2011), in which day of capture (1 = 1 Jan 2017) was included as the smoothed term predicting the probability of infection. These models assess the finer variation in infection throughout the season for the various parasites, but do not allow us to test interactions between smoothed terms to directly determine whether parasites have different seasonalities. In order to compare the seasonality of infection probability between lineages, species and PCR methods, we used GLMs with logit-link function and binomial distribution, but the response variable was defined as the number of infected/number of non-infected individuals. In these analyses, we specifically tested for an interaction between sampling period and lineage/species/method, with country included as a predictor whenever appropriate. Significance of the variables was assessed using Chi-square tests with the ANOVA function.

Data accessibility

Raw data are available in Mendeley Data, DOI: https://doi.org/10.17632/hf8kksbjv8.1.

Results

Parasite prevalence and diversity varies between sites

A total of 21 haemosporidian lineages were found in 733 blood samples (including retraps) of house sparrows collected in Spain, Bulgaria, Poland and Sweden: 11 Plasmodium, seven Leucocytozoon, of which two were new (L_PADOM36 and L_PADOM37, GenBank accession numbers MT080104 and MT010105, respectively), and three Haemoproteus. The number of lineages per country declined slightly towards the north, and the parasite communities at each site differed (Table 1, Table 2, Supplementary Table S1). Only the most common lineage (P_SGS1, belonging to P. relictum; see Fig. 1) occurred at all four study sites. Furthermore, the second most common lineage (H_PADOM05, belonging to H. passeris; Fig. 1) occurred almost exclusively in Spain; whereas Leucocytozoon was only detected in Sweden and Spain, and these sites did not share any lineage of this genus (Table 1, Table 2).

Table 1

Number of lineages and overall prevalence of haemosporidian parasites in house sparrows, detected using the nested PCR method.

	Sweden	Poland	Bulgaria	Spain	Total
Overall prevalence (%)	13.27	37.72	50.58	79.29	48.02
Sample size (incl. retraps)	196	167	172	198	733
Total number of lineages	6	8	8	9	21
Plasmodium	3	6	7	3	11
Leucocytozoon	3	0	0	4	7
Haemoproteus	0	2	1	2	3

Table 2

Number of house sparrows (including retraps) for which each lineage and lineage combination of haemosporidian parasite was detected using the nested PCR method. See also Supplementary Table S1 for the number of samples per lineage and country.

		Sweden	Poland	Bulgaria	Spain	Total
	Uninfected	170	104	85	41	381
	Unidentified P/H	0	2	0	1	3
Plasmodium	SGS1	9	35	56	14	114
	SGS1/COLL1		2	4		6
	SGS1/GRW11			1		1
	SGS1/PADOM01			1		1
	BT7	2				2
	COLL1		3	3	1	7
	GRW06			1		1
	GRW11		9	12		21
	GRW11/PADOM02			2		2
	GRW11/COLL1			1		1
	LINOLI01		1			1
	PADOM02		3	1		4
	PADOM2/SYAT24			1		1
	PAGRI02				1	1
	SYAT24			2		2
	TURDUS1	3	1			4
P/L	SGS1/RECOB3				2	2
	SGS1/SYCON06				1	1
Leucocytozoon	HIRUS07				1	1
	BT2	5				5
	PARUS21	5				5
	PADOM36/BT2	1				1
	PADOM36	1				1
	RECOB3/SYCON06				1	1
	RECOB3				6	6
H/L	PADOM22/RECOB3				1	1
	PADOM05/PADOM37				1	1
	PADOM05/SYCON06				1	1
	PADOM05/HIRUS07				1	1
	PADOM05/RECOB3				11	11
	PADOM5/PADOM22/RECOB3				1	1
H/P	PADOM05/PAGRI02				1	1
	PADOM05/PADOM02		1			1
	PADOM03/SGS1		1			1
Haemoproteus	PADOM05			2	100	102
	PADOM05/PADOM22				11	11
	PADOM22				1	1
	PADOM03		5			5
	Total	196	167	172	198	733

The total prevalence differed markedly between the four countries, from 79% in Spain to 13% in Sweden (Table 1). Our initial screening of haemosporidian parasites, using only the nested PCR method, resulted in low numbers of co-infections between Plasmodium and Haemoproteus lineages (only three) despite the high number of single infections of each genus (Table 2). Moreover, the number of Plasmodium lineages detected in Spain was low considering the expected general increase in diversity with decreasing latitude (Table 1). We therefore proceeded to investigate whether the high prevalence of H_PADOM05 in Spain concealed the presence of Plasmodium parasites.

Nested PCR underestimates the occurrence of co-infections

The use of a multiplex PCR that reports the occurrence of three genera (Haemoproteus, Plasmodium and Leucytozoon) was compared with the prevalence of infection resulting from the nested PCR in 215 samples shown to be infected by the latter method. Only 60% of the samples agreed fully, hence these screening methods produced inconsistent, and biased, results for many of the tested samples (Table 3). Indeed, there was a high number of samples (75 out of 131; 57.2%) for which the nested PCR only detected Haemoproteus lineages (or Haemoproteus and Leucocytozoon), but the multiplex PCR detected the presence of both Haemoproteus and Plasmodium lineages (Table 3). This contrasts with the total agreement between the methods when the nested PCR detected only Plasmodium or Plasmodium and Leucocytozoon (n = 63, see Table 3). In addition, there were only 10 out of 78 cases (12.8%) in which the multiplex PCR method detected Plasmodium (or Plasmodium and Leucocytozoon) and the nested PCR method detected Haemoproteus. That is, the failure of the nested PCR to detect Plasmodium in the presence of Haemoproteus is significantly higher than both: (i) its failure to detect Haemoproteus in the presence of Plasmodium (Chisq = 56.4; d.f. = 1; P < 0.001), and (ii) the failure of the multiplex PCR to detect Haemoproteus in the presence of Plasmodium (Chisq = 38.2; d.f. = 1; P < 0.001).

Table 3

Comparison between nested and multiplex PCRs regarding the detection of parasites belonging to the genera Haemoproteus (H), Plasmodium (P) and Leucocytozoon (L) in single and co-infections. Numbers represent the number of samples, the diagonal represent the cases when the methods were consistent, and those in red/bold are the number of cases when Plasmodium was not detected by the nested PCR but was present according to the multiplex PCR in co-infections with Haemoproteus.

In order to determine the species of the hidden Plasmodium infections, we designed Plasmodium-specific primers and did a third molecular genetic screening. In 43 out of 67 samples for which the nested PCR only detected Haemoproteus (or Haemoproteus and Leucocytozoon) the newly-designed, Plasmodium-specific primers confirmed the presence of lineages belonging to the P. relictum clade, the P. cathemerium clade or, in nine cases, mixed infections of lineages belonging to both of these clades (see Fig. 1). Notably, two of the nine samples that had mixed infections of Plasmodium lineages also had H_PADOM05 and L_RECOB3. That is, these two birds were simultaneously infected by at least four haemosporidian lineages. It is also worth noting that one sample that had P_SGS1 according to the nested PCR method and only Plasmodium according to the multiplex PCR method, actually had a mixed infection of P. relictum and P. cathemerium, as revealed by the newly-designed primers. Hence, it is clear that co-infections are underestimated by the nested PCR method, especially those including Haemoproteus and Plasmodium.

Importantly, the patterns of seasonal occurrence of Plasmodium in Spain depended on the molecular method used. According to GAMs, when only the nested PCR was considered, the probability of being infected by Plasmodium parasites in Spain (i.e. mostly P. relictum) appeared relatively high during the winter period, but then decreased markedly towards summer (as Plasmodium infections were missed by the nested PCR when H_PADOM05 became very abundant in the summer months; see Fig. 2). However, when the results from the multiplex PCR were also taken into consideration, the pattern of occurrence of Plasmodium parasites was much more flat (Fig. 2). This difference in seasonality is supported by a significant interaction between sampling period and PCR method affecting the probability of infection (GLM Chisq: PCR Method: P < 0.001, Sampling Period: P = 0.009; PCR Method versus Sampling Period: P < 0.001; see Supplementary Fig. S1).

Fig. 2

Seasonal variation of probability of infection by Haemoproteus and Plasmodium spp. in Spain. Curves depict coefficients and confidence intervals of three separate Generalized Additive Models, one of which analyses the results obtained for Plasmodium using nested PCR only (see also Supplementary Fig. S1). For statistics see Section 3.

Seasonality varies between lineages and sites

No significant differences in overall infection probability (i.e. for all parasites combined) were detected between the sexes, controlling for country and sampling period (GLM Chisq test: P = 0.92), even when Spain (where the overall probability of infection was very high) and Sweden (where probability of infection was low) were excluded (all P > 0.05). Similarly, there were no differences in infection probability between age classes in July and September (i.e. when age can be reliably determined; GLM Chisq test: P = 0.50). Consequently, age and sex were not included in the analyses of seasonal variation.

Overall infection probability (i.e. for all parasites combined) differed significantly between countries, and there was a highly significant interaction between country and sampling period, indicating that seasonality in prevalence is site-dependent (GLM Chisq tests: Country: P < 0.001; Sampling Period: P = 0.059; Country * Sampling Period: P < 0.001). This is because prevalence was higher in Spain, decreasing towards the north (Table 1), and seasonality was also more pronounced in Spain (see Fig. 3). Separate analyses for each country showed that overall infection probability varied highly significantly with sampling period in Spain (GLM Chisq test: P < 0.001) and tended to do so in Poland (P = 0.067). However, the overall infection probability did not differ between sampling periods in either Bulgaria (P = 0.271) or Sweden (P = 0.109).

Fig. 3

Seasonal and geographical variation of haemosporidian parasites in house sparrows. Results are from the nested PCR for first captures only (i.e. retraps were excluded). No sampling took place in Bulgaria during September 2017.

Analyses of the most abundant lineages and species showed a variety of patterns in seasonal variation (Fig. 4). H_PADOM05 showed a strong seasonality in Spain (GLM Chisq test P < 0.001; GAM Chisq test P < 0.001), where the infection probability was low in the winter months and then increased progressively until reaching a peak in autumn (Fig. 2). In contrast, the probability of being infected by the commonest Leucocytozoon lineage (L_RECOB3), which occurred exclusively in Spain (Table 2, Supplementary Table S1), did not vary significantly with time (GLM chisq P = 0.278; GAM chisq = 0.90; P = 0.645; Fig. 3D). On the other hand, the probability of infection by P. relictum (i.e. P_SGS1 and P_GRW11; see Fig. 1) varied between countries and sampling periods (GLM Chisq test: Country: P < 0.001; Sampling Period: P = 0.004) and, interestingly, the pattern of seasonal variation varied between countries (Country versus Sampling Period: P = 0.002). The same pattern occurred when the analyses were restricted to P_SGS1 (GLM Chisq test: Country: P < 0.001; Sampling Period: P = 0.06; Country versus Sampling Period: P = 0.009; Spain being excluded as P_SGS1 could not be positively identified from co-infections). Differences in seasonality between countries is supported by GAMs showing that the probability of infection by P. relictum varied significantly with day of the season in Spain (GAM Chisq = 13.45, P = 0.017) and in Sweden (GAM Chisq = 3.91, P = 0.048), but not in Poland (GAM Chisq = 10.96, P = 0.109) nor Bulgaria (GAM Chisq = 8.61, P = 0.128). The infection probability of P. relictum generally declined from September 2016 to September 2017 (Fig. 4A), showing a slight increase in April only in Spain (Fig. 4B). In contrast, the probability of infection by P. cathemerium (lineages P_COLL1, P_PADOM01 and P_PADOM02; see Fig. 1) did not seem to differ between countries, but approached significance for sampling period and the interaction of these two variables despite its smaller sample size (GLM Chisq test: Country: P < 0.001; Sampling Period: P = 0.06; Country versus Sampling Period: P = 0.07). A GAM for P. cathemerium with country included as a predictor showed a significant effect of day of capture on the probability of infection (GAM Chisq = 12.26, P = 0.031), and a pattern of seasonal variation with two peaks of occurrence in September (higher in 2016 than 2017) and April/May (Fig. 4C).

Fig. 4

Seasonal variation of infection by Plasmodium relictum (P_SGS1 and P_GRW11) for (A) all countries combined and for (B) Spain only, by (C) Plasmodium cathemerium (P_COLL1, P_PADOM01 and P_PADOM02) for all countries and by (D) L_RECOB3, which was found only in Spain, from September 2016 until September 2017 (sampling occurred every second month). These results depict the coefficient of probability of infection and confidence intervals of Generalized Additive Models in which day of capture was included as the smooth variable. For statistics see Section 3.

The probability of infection by the two lineages within P. relictum (P_SGS1 and P_GRW11) differed significantly from each other and between countries (with Spain excluded as these lineages could not be identified in mixed infections with the newly-designed primers), and almost significantly between sampling periods, but seasonality differences between countries were not detected (GLM Chisq test: Lineage: P < 0.001; Country: P < 0.001; Sampling Period: P = 0.06; Lineage versus Sampling Period: P = 0.885). Comparisons of the probability of infection between P. relictum and P. cathemerium did not produce significant interactions between sampling period and species (GLM Chisq test: Species: P < 0.001; Country: P < 0.001; Sampling Period: P = 0.003; Species versus Sampling Period: P = 0.385), indicating that the pattern of seasonal variation by these taxa was similar when controlling for differences between countries. However, the probability of infection by Haemoproteus (mostly H_PADOM05) and Plasmodium (mostly P. relictum) interacted significantly with season in Spain (GLM Chisq test: Genera: P = 0.019; Sampling Period: P < 0.001; Genera versus Sampling Period: P < 0.001). In addition, there was a highly significant interaction between sampling period and parasite when comparing Haemoproteus and L_RECOB3 infections is Spain (P < 0.001); and an almost significant interaction between the same variables when comparing P. relictum and L_RECOB3 (P = 0.054); indicating that all these three lineages vary in their seasonal patterns of occurrence in Spain.

Discussion

Our study on seasonal variation of haemosporidian parasites in resident house sparrows benefited from concurrent sampling at multiple, distant sites within the European temperate region. This allowed us to compare the prevalence and diversity of these parasites between locations, as well as the patterns of seasonal variation between the most common lineages, species and sites. This setup enabled us to avoid the bias that can be caused by common confounding factors such as sampling year, the population origin of hosts and grouping of parasite lineages that are distantly related. In addition, the use of multiple genetic methods allowed us to identify and overcome limitations that potentially could lead to incorrect ecological inferences.

We observed the expected decline in prevalence and diversity of haemosporidian parasites with increasing latitude (Table 1). The observed diversity of parasites, however, is strongly dependent on the sample size of infected hosts, particularly in haemosporidians for which there can be a very large number of lineages infecting the bird community at a single site and spill-over appears common (Neto et al., 2015, Ellis et al., 2020). Although the number of house sparrows sampled in each country was relatively similar, the number of infected birds (i.e. those that can actually be used to determine parasite diversity) was very different, as the overall prevalence of malaria varied from 79% in Spain to only 13% in Sweden (Table 1). Furthermore, as outlined in the results, it is likely that some Plasmodium lineages were missed in our Spanish sample due to the underestimation of mixed infections with Haemoproteus parasites (see below), and the inability of our Plasmodium-specific primers to identify all lineages as defined in the MalAvi database (Bensch et al., 2009). Hence, the latitudinal comparison of the diversity of parasites should be considered provisional, and additional sites should eventually be included to randomize site-specific effects (e.g. proximity to water and temperature; Wood et al., 2007; Lachish et al., 2011; Loiseau et al., 2013; Ferraguti et al., 2018). Importantly, given the differences in overall prevalence and in the actual parasite composition between sites in our study, the selection pressure exerted by the parasites on house sparrows’ immune system could conceivably lead to the evolution of local adaptation, which is the focus of our forthcoming studies.

The unexpected low diversity of Plasmodium lineages in Spain, together with the high prevalence of a Haemoproteus lineage, led us to believe that the classical nested PCR method widely used in studies of avian haemosporidian parasites may have failed to detect mixed infections, which has been reported previously in the literature (Valkiūnas et al., 2006, Bernotienė et al., 2016, Ciloglu et al., 2019). Our comparison between the nested PCR and the multiplex PCR revealed multiple cases of inconsistencies that were biased towards failure to detect Plasmodium in the presence of Haemoproteus parasites by the nested PCR. The reason for this seems to be that the parasitaemia of Haemoproteus is generally much higher than the parasitaemia of the Plasmodium parasites (as observed in most of our samples that were tested with quantitative PCR; unpublished results), making it very difficult to visually detect ambiguities in the electropherograms. The development and use of Plasmodium-specific primers further confirmed that Plasmodium lineages had been missed in many samples with H_PADOM05 and other Haemoproteus parasites, and made it possible to detect individual birds that had simultaneous infections of up to four haemosporidian parasite lineages in their blood. We find it likely that multiple infections of parasites are the norm rather than the exception, at least in some hosts and at some locations, especially if multiple tissues are analysed (Svensson-Coelho et al., 2016). Although there are additional protocols and primers that were not tested in our study (e.g. Bernotienė et al., 2016, Pacheco et al., 2018), current molecular methods are limited in the number of detectable lineages and/or are too time consuming (e.g. cloning). Hence, the use and development of new methods (e.g. metabarcoding or even metatranscriptomics; see Galen et al., 2020) could potentially be used in the future, for instance, to study the accumulation of parasites over a bird’s lifetime, facilitation/competition between lineages, as well as to obtain more accurate measures of prevalence of the various lineages. It is also worth pointing out that positive PCR detection of parasites in the blood does not mean that these parasites are available for transmission: microscopy can be important to complement molecular studies, for instance, to confirm the presence of gametocytes.

Our results further show that at least in some circumstances, particularly when Haemoproteus infections are frequent, the sole use of the nested PCR method can lead to false ecological inferences. This is because the seasonal variation of Plasmodium parasites in Spain significantly depended on the molecular genetic method used (Supplementary Fig. S1). When only the nested PCR results were considered, it appeared as if the probability of being infected by Plasmodium was higher during winter and declined during spring/summer, when the Haemoproteus infections increased (Fig. 2, Supplementary Fig. S1). This could potentially be interpreted as a case of competition between the parasites. However, the inclusion of results from the multiplex PCR produced a significantly different, more flat pattern of variation of the Plasmodium parasites during the season in Spain (Fig. 2, Supplementary Fig. S1). Therefore, we recommend that special attention should be paid to such biases, in which case multiple molecular methods should be used.

Our analyses revealed the existence of both site-specific and lineage-specific patterns of seasonal variation of haemosporidian parasites in house sparrows. Cosgrove et al. (2008) had already shown that Plasmodium spp. can differ in seasonality within the same site, so different patterns of seasonal variation between sites are not surprising, especially if the parasite communities vary between sites, as in this case. Surprisingly, the haemosporidian parasites showed greater differences in prevalence across the season in Spain than in the more northern countries, although seasonality is expected to increase with latitude (Cook et al., 1990, Dowell, 2001, Lisovski et al., 2017; but see Altizer et al., 2004 for an exception). However, this pattern resulted mostly from the strong seasonality of H_PADOM05, which occurred almost exclusively in Spain, but partly also from differences in seasonality of P. relictum between the countries. Plasmodium relictum seemed to show a double peak of occurrence in Spain, but in the other three countries there was a flat or a general decline in prevalence throughout the sampling period (Fig. 4). This decline seems to have resulted mostly from differences between years rather than from seasonal variation, as the prevalence (of both P. relictum and P. cathemerium) was much higher in September 2016 than in September 2017 (which may be associated with differences in weather and availability of vectors between years; see Diouf et al., 2017). Seasonality of P. relictum is also known to vary along the altitudinal gradient in the Hawaiian Islands, where there is a high transmission and minor seasonal or annual variation in low-elevation forests (as was also described by Cosgrove et al., 2008 in the U.K.), episodic transmission in mid-elevation forests with site to site, seasonal and annual variation, and only a slight risk of infection during summer in high elevation forests (Samuel et al., 2011).

Comparisons between lineages (P_SGS1 and P_GRW11) within the same species (P. relictum) were constrained by a much lower prevalence of GRW11, as were comparisons between P. relictum and P. cathemerium, the latter having much lower prevalence. It is possible that these lineages and species actually differ in probability of infection throughout the season, but we were unable to detect these differences in GLMs due to low statistical power and/or relatively rough temporal resolution. However, GAMs indicate that P. cathemerium appeared to have a double peak of occurrence when all countries were combined, the first in April/May and the second during autumn, whereas this double peak pattern was only apparent for P. relictum in Spain. Nevertheless, we detected clear differences in seasonal variation between P. relictum/P. cathemerium, H_PADOM05 and L_RECOB3. Indeed, H_PADOM05 showed a very strong single peak of prevalence in the early autumn, whereas L_RECOB3 showed a flat pattern of occurrence.

The double peak in prevalence of P. relictum in Spain and for P. cathemerium for all countries combined, and the annual differences in their autumn prevalences, seem to match the model described by Beaudoin et al. (1971). However, because these species infect a huge diversity of bird hosts (see MalAvi database), it is very difficult to infer the causes of variation in prevalence during the season. Mosquitos are already on the wing in April/May in Spain, when the first peak in prevalence occurs, and there are also plenty of migratory birds of other species that have already arrived from the south, potentially being the sources of parasites for new transmissions. It is not clear whether the spring peak occurred due to relapses or partly also due to new infections, although the autumn peak included many juvenile birds and thus surely involved new transmission to susceptible individuals. The other two abundant parasites (H_PADOM05 and L_RECOB3) did not show the same pattern of occurrence as described by the model in Beaudoin et al. (1971), which was in fact described for Plasmodium only. Nevertheless, the same processes (relapse, transmission to juvenile hosts) may be operating, but because the vectors of these lineages are not known, it is not possible to determine at which time of year the parasites should invest in transmission (Pigeault et al., 2018), nor whether this coincides in time with the production of juvenile birds (in the case of H_PADOM05). Culicoides biting midges are known to transmit Haemoproteus lineages, and at least some species of these insects seem to be generalists relative to the parasites (Valkiūnas, 2005, Chagas et al., 2019). However, this might not always be the case, and many species of Culicoides vary dramatically in their patterns of seasonal and geographic abundance (Cuéllar et al., 2018), making predictions difficult. Clearly, more effort should be put into describing the relationships between haemosporidian parasites and their vectors (e.g. Ferraguti et al., 2013; Chagas et al., 2019).

In conclusion, our data is well in line with observations of different seasonalities between different parasite lineages/species by Cosgrove et al. (2008), but we also find patterns of seasonal occurrences that agree with the one described by Beaudoin et al. (1971). Importantly, we also show the existence of different seasonalities of the same parasite lineage between countries. It is not possible to infer the causes of variation in seasonal patterns of occurrence, as very little is known about: (i) the vectors transmitting many of the lineages, (ii) the probability of clearing infections of particular parasites by the hosts, nor (iii) the role of co-occurring bird species as reservoirs of some parasites. It is nevertheless evident from our results that pooling different parasite lineages (even within the same genus), years and sites in statistical analyses is likely to obscure their patterns of occurrence (which are differentially influenced by the relative abundance of the parasites), making ecological inferences difficult. Pooling parasites, in particular, should preferably involve testing differences between lineages/species, and acknowledgment that the patterns are mostly driven by the commonest lineages, avoiding extrapolation of the results to the whole species or genus. Finally, limitations of the molecular screening methods should be considered and evaluated, and when biases are found the use of multiple molecular methods is recommended.