Coevolution between hosts and parasites with partially overlapping geographic ranges

Introduction

The majority of interacting species have geographical ranges that do not completely overlap. Many host and prey species, for instance, interact with a particular parasite or predator within only a fraction of their geographical ranges. Examples arise within a broad range of natural, agricultural, and medically important interactions, including malaria and humans (Ruwende et al., 1995), lodgepole pine and pine squirrels (Benkman, 1999), magpies and cuckoos (Soler et al., 1999), salamanders and predatory fish (Storfer & Sih, 1998), powdery mildew and wheat (Paillard et al., 2000), and snails and trematodes (Lively & Jokela, 1996). In addition, interactions between partially overlapping species readily arise as a consequence of biological invasions. The result is often a mosaic of coevolutionary hot spots, where host and parasite overlap and have reciprocal fitness effects, and cold spots, where the parasite is absent.

The geographic mosaic theory of coevolution suggests that variation in selection across coevolutionary hot and cold spots is of crucial importance to evolving interactions. This conceptual framework argues that much of the dynamics of coevolution result from spatial variation between coevolutionary hot and cold spots, in conjunction with selection mosaics, gene flow, genetic drift, and metapopulation dynamics (Thompson, 1994). These ideas have been incorporated into recent mathematical models investigating how the interplay between these factors shapes evolving interactions. This work has included studies of gene-for-gene or matching alleles coevolution within metapopulations, (e.g. Gandon et al., 1996; Gandon, 1998, 2002; Burdon & Thrall, 1999; Damgaard, 1999; Parker, 1999; Thrall & Burdon, 1999, 2002), variable patterns of coevolutionary selection (Nuismer et al., 1999, 2000; Gomulkiewicz et al., 2000), and models that incorporate spatial variation in demographical parameters (Hochberg & van Baalen, 1998; Hochberg et al., 2000). Together, these models have suggested that selection mosaics, gene flow, and metapopulation dynamics may play an important role in evolving interactions between completely overlapping species.

The results of several of these studies have also indicated that geographically structured coevolution could be important in generating observed patterns of local adaptation in host–parasite interactions. The work of Gandon et al. (1996), for instance, used simulations to demonstrate that spatial patterns of local adaptation often emerge as a result of differences in the phase of allele frequency oscillations in different populations. In addition, this work demonstrated that the extent of parasite local adaptation depends upon the relative rates of gene flow in host and parasite. Similar results were found by Lively (1999), who analysed an extreme scenario where only the parasite moved, and by Gandon (2002) through analytical approximations. Using a different theoretical framework, Nuismer et al. (2000) showed that asynchronous cycles and local adaptation were often ephemeral in the presence of gene flow. Only when there was substantial spatial heterogeneity in parasite virulence could asynchronous cycles be maintained. Together, these studies have demonstrated that the rate of gene flow is a crucial determinant of emerging patterns of local adaptation in parasites.

An assumption common to all of these models, however, is that the geographical ranges of both species coincide and, with the exception of transient local extinctions, coevolve everywhere. Although these assumptions are surely met in many interactions, recent empirical studies have demonstrated that other interactions may be better characterized as a stable mosaic of coevolutionary hot and cold spots (Thompson, 1994; Brodie & Brodie, 1999; Benkman et al., 2001; Brodie & Ridenhour, 2002; Parchman & Benkman, 2002; Thompson & Cunningham, 2002). This facet of the geographic mosaic theory has received significantly less attention than those cases where both host and parasite coevolve and coexist throughout their geographical ranges. At this point only one model has explored spatial variation between coevolutionary hot and cold spots (Gomulkiewicz et al., 2000). This model analysed a case where both species overlap completely, but coevolve within only a subset of the available habitat. Despite its relative simplicity, the results of this model demonstrated that spatial variation among coevolutionary hot and cold spots could produce novel coevolutionary dynamics.

Here we extend the theory of coevolutionary hot spots to incorporate explicit spatial structure and interactions between a host and parasite with partially overlapping geographical ranges. Our approach is to develop a spatially explicit genetic model of coevolution between a host and a partially overlapping parasite. To the extent that our model couples coevolutionary hot and cold spots with gene flow, the model is similar to that of Gomulkiewicz et al. (2000). Our model, however, is quite different in its capacity to address explicit spatial patterns of allele frequencies and maladaptation. This feature of the model allows us to explore the following unanswered questions. How does the explicit spatial location of coevolutionary hot spots modulate their impact on evolving interactions? Do foci of maladaptation emerge within the geographical range of an interaction, and if so, where? What conditions must be satisfied for coevolutionary cycles to propagate? In addition to answering these specific questions, our results directly extend previous models in an effort to continue development of a mathematical framework for the geographic mosaic theory of coevolution.

The model

We assume that the host is distributed among a linear series of populations. In the simplest possible spatial arrangement, a coevolutionary hot spot of length H abuts a cold spot of length C (Fig. 1a). Populations of the parasite are assumed to exist within a contiguous subset of these localities labelled −H to 0. Within the region of species overlap (−H to 0), selection imposed on the two species by the interaction is reciprocal. By definition, this geographical region constitutes a coevolutionary hot spot (Thompson, 1994; Gomulkiewicz et al., 2000). Outside of this region, in localities labelled 0 to C, the parasite is absent. In addition to this simple configuration, we will consider a scenario where a coevolutionary hot spot is located between two cold spots (Fig. 1b).

We assume that both species are haploid, with alleles Y and y in the parasite and Z and z in the host. Within hot spots, fitnesses of both species are determined by a matching alleles model such that encounters between Y and Z or y and z result in a fitness decrement to the host and a fitness increment to the parasite (e.g. Seger, 1988; Gavrilets & Hastings, 1998). In both hot and cold spots the host is assumed to experience directional selection favouring the z allele, perhaps because of an associated cost or pleiotropic effect as is the case for some forms of malaria resistance (Ruwende et al., 1995). We model these assumptions with the following fitnesses W_k,i for each genotype k, at geographical location i:

((1a))

((1b))

((1c))

((1d))

In these expressions p_1,i is the frequency of the Y allele in parasites and p_2,i is the frequency of the Z allele in hosts, s is the coefficient of directional selection imposed on hosts by factors extrinsic to the interaction, c is the fitness sensitivity of the parasite to changes in host allele frequencies, and b is the fitness sensitivity of hosts to changes in parasite allele frequencies within hot spots. Without loss of generality, we assume that s,b,c > 0. These fitness functions, when combined with the standard equations of population genetics result in a system of two recursion equations describing allele frequencies in both species after one round of selection (Seger, 1988; Gavrilets & Hastings, 1998; Nuismer et al., 2000).

We incorporate gene flow as in Nuismer et al. (2000), by assuming that after selection each species moves between habitats with a rate determined from a discretized zero mean Gaussian distribution. Thus the probability of movement between patches decreases as the distance between patches increases. The rates of movement in the two species can differ and are determined by the parameters for the parasite, and for the host. These parameters measure the dispersal variances of the two species such that increasing values correspond to increasing movement between patches. Making the additional assumption that selection is weak, results in the following system of approximate partial differential equations:

((2a))

((2b))

A more detailed development of these reaction–diffusion approximations, as well as discrete analogues are found in Nuismer et al. (2000). Here p₁(x,t) and p₂(x,t) are the allele frequencies in species 1 and 2 at time t, and the variable x is the continuous analogue of the position variable i in the discrete model (Nuismer et al., 2000).

We compared analytical solutions of equations (2) to numerical simulation of the discrete model, and found them to be in good agreement whenever selection is sufficiently weak (e.g. s,b,c < 0.05). For equations (2), we assume zero flux boundary conditions, such that no genes cross the left or right hand boundaries, and that, for hosts, allele frequencies and their slopes are continuous at x = 0. As long as spatial variation in selection occurs along only one spatial dimension, the model is equally applicable to a two dimensional space.

Parasite at the periphery of host range

Using the model described above, we evaluated how coevolution between a host and a peripherally overlapping parasite influences the overall evolution of the interaction. As a null expectation, we first considered the case where the parasite is completely absent (H = 0). Under these conditions, directional selection leads to fixation of the z allele in the host, resulting in the stable monomorphic equilibrium

(3)

If, however, the host overlaps and coevolves with the parasite within a fraction of its geographical range (H > 0), fixation of this allele is no longer inevitable.

Consider the equilibrium where the host is fixed for the z allele favoured by directional selection, and the parasite is fixed for the matching y allele:

((4a))

((4b))

At this ‘global matching’ equilibrium the parasite is perfectly adapted to the host while the host is perfectly maladapted to the parasite. For a coevolutionary hot spot to be of consequence, it must destabilize this global matching equilibrium.

We determined the conditions under which coevolution destabilizes this equilibrium by performing a local stability analysis using the continuous approximation (2) following the methods used in Nuismer et al. (2000). This analysis shows that equilibrium (4) is unstable if

(5)

where θ₀ is the unique root between 0 and π/2 of

(6)

with γ = C/H, and α = |s/(b−s)|. This result has two important biological consequences. First, it demonstrates that, unless the strength of parasite selection on the host exceeds the strength of directional selection (i.e. b > s), no amount of overlap between species will lead to coevolution. Secondly, it shows that even when this is true, a critical threshold of overlap must still be exceeded for coevolution to destabilize the global matching equilibrium favoured by directional selection. Specifically, this critical threshold decreases as rates of host gene flow decrease (Fig. 2a, left hand line) and the strength of coevolutionary selection increases.

When both species overlap sufficiently for (5) to be satisfied, simulation of the discrete model reveals that static polymorphic clines evolve within the host and potentially the parasite species as well (Fig. 3). Two features of these clines stand out. First, polymorphism is focused within the region where both species overlap. Secondly, at least for relatively large regions of overlap, these clines are nonmonotonic (Fig. 3b). This latter feature is a common signature of coevolutionary clines (Nuismer et al., 2000) not generally found in single species models with gene flow and spatially variable selection (Slatkin, 1973; Nagylaki, 1975; Garcia-Ramos & Kirkpatrick, 1997).

We analysed how clines affect the spatial structure of maladaptation within host and parasite by defining an index of maladaptation for species j in population i as the proportion

(7)

where , and , are the maximum, minimum, and observed population mean fitness of species j in population i respectively. Thus populations with the minimum possible local mean fitness have a maladaptation index of 1, whereas populations with the maximum possible local mean fitness have a maladaptation index of 0. This measure of maladaptation is quite different from the measure of local maladaptation often used in studies of host–parasite coevolution (e.g. Gandon et al., 1996; Kaltz et al., 1999; Gandon & Michalakis, 2002). Instead of defining maladaptation relative to other populations (e.g. Gandon, 2002), (7) defines maladaptation in terms of the theoretical maximum and minimum values that population mean fitness can take at a given locality. Thus (7) is a measure of the local ‘coevolutionary load’ a population experiences, and is independent of the genetic state of other populations. In contrast, local maladaptation depends upon the genetic state of other populations and can only exist if gene frequencies are spatially heterogeneous (Nuismer et al., 2000; Gandon, 2002). In what follows, we will refer to maladaptation relative to other populations as local maladaptation, and maladaptation given by (7) as maladaptation.

Plotting the value of maladaptation determined from (7) across space reveals that host maladaptation is concentrated at the periphery of the parasites geographical range (Fig. 3, right hand panels). In contrast, parasite maladaptation is minimized within this same region (Fig. 3, right hand panels). Simulations show these trends to be very general, occurring across many combinations of parameter values. At the same time, the results of simulations suggest that substantial local maladaptation is present in both species as well. This local maladaptation arises as a consequence of the spatial genetic structure that is maintained in the form of static allele frequency clines (Fig. 3, left hand panels).

We further investigated the dynamics of maladaptation using simulations to calculate average levels of maladaptation across space in the parasite and host for different values of gene flow. Average host maladaptation increases with host gene flow, and the fraction of habitat occupied by the parasite (Fig. 4, circles). This result is quite different from that of Gandon (2002), who found that increased host gene flow often caused improved host local adaptation to evolve. The reason for this discrepancy is that here, maladaptation evolves simply because hosts must simultaneously adapt to both the parasite and the abiotic environment in the face of gene flow. Much more interesting, however, are average levels of maladaptation confined to the geographical region where both species overlap. Within this region, host and parasite maladaptation depends upon the fraction of the host's geographical range occupied by the parasite. In particular, the parasite species becomes increasingly maladapted as it occupies a larger fraction of its host's geographical range (Fig. 4, squares). In contrast, as the parasite occupies an increasing fraction of the host's geographical range, the host becomes increasingly well adapted to the parasite (Fig. 4, triangles). These transitions in the level of reciprocal adaptation occur quite rapidly. Thus to the extent that a parasite's geographical range is limited by adaptation to its host, this mechanism may play a role in limiting the geographical range of a parasite.

If the parasite species is present throughout much of the host's geographical range, or if the selective impact of parasites on hosts is large, permanent oscillations in allele frequency may evolve. Numerical simulation of the discrete model suggests that two conditions must be met. First, the overall strength of parasite selection on the host must exceed that of directional selection:

(8)

Secondly, host gene flow must exceed some critical threshold value, with this value decreasing with the fraction of host range occupied by the parasite (Fig. 2a, right hand line) and with decreases in the strength of directional selection in cold spots.

A striking result from the simulations is that under these conditions allele frequency oscillations are perfectly in phase across the range of the interaction. Under these conditions, local adaptation of host and parasite is absent as both species are genetically homogenous across space. For the values of selection we considered (e.g. s,b,c < 0.1), we were unable to find parameter values consistent with the permanent and asynchronous oscillations (transient asynchronicity, however, can persist for many thousands of generations) required for dynamic patterns of local adaptation (e.g. Gandon et al., 1996; Morand et al., 1996). Consequently, if selection is not extremely strong, permanent patterns of local adaptation must be associated with static allele frequency clines rather than asynchronous oscillations. This result provides a direct prediction for studies of the dynamics of local adaptation in natural populations.

Parasite at the centre of host range

In this section we consider how the explicit spatial location of overlap between species modulates the impact of a coevolutionary hot spot. We repeated several of the analyses from the previous section, but assumed the parasite species occurs within the centre of the host's geographical range (Fig. 1b). Because we assume the habitat is symmetrical about the region where species overlap occurs, we could analyse this case by simply reflecting the original axis about the origin at x = 0 and re-scaling the parameters H and C such that the sizes of the hot and cold spot habitats remain constant (Fig. 1b) (e.g. Slatkin, 1973). Provided the symmetry of this re-scaling is maintained, the results derived in this section also apply to a periodic environment where coevolutionary hot and cold spots alternate over space (For the single species case see Slatkin, 1973; Roughgarden, 1996). For all cases, we compared these new results to those from the previous section, where the species overlapped only at one edge of the host's geographical range.

We first derived conditions that allow parasite driven selection on the host to destabilize the equilibria favoured by directional selection. Repeating the stability analysis of the first section, except with the spatial coordinates depicted in Fig. 1b, shows that the global matching equilibrium favoured by directional selection (4) is unstable only if

(9)

where θ₀ is defined by (6). Comparing (5) and (9) for any fixed set of gene flow and selection parameters demonstrates that a centrally located parasite species must be present in a larger fraction of the host's geographical range for global matching to be destabilized. Biologically, this means that a greater amount of overlap between species is required to elicit a coevolutionary response in the host when overlap occurs at the centre vs. the periphery of the host's range (Fig. 2b, left hand line). In contrast, simulations suggest that when condition (8) is satisfied and allele frequency oscillations become possible, oscillations propagate more readily when the parasite is located centrally. This effect can be most readily visualized by comparing the right hand lines of Fig. 2a and b: for any value of host gene flow, the critical fraction of overlap required for permanent oscillations is substantially reduced when the parasite is centrally located.

Together, these results make good intuitive sense. When a parasite is confined to the centre of the host's geographical range, host populations encountering parasites receive gene flow from populations where the parasite is absent on both sides. Alternatively, if the species overlap only at the edge of the host's geographical range, host populations encountering parasites receive only half this amount of gene flow from populations on the side where the parasite is absent. Consequently, when parasites exert selection on hosts that is weak relative to directional selection [bH < s(C + H)], the coevolution of clines is considerably more likely when species overlap occurs in peripheral, rather than central, host populations. When selection exerted on the host by the parasite is strong relative to directional selection [bH > s(C + H)], however, global allele frequency oscillations are more likely when coevolutionary hot spots are located within the centre of the host's geographical range. This is a result of the greater efficacy with which a centrally located coevolutionary hot spot can project migrants into adjacent cold spots.

Discussion

Our analyses show that coevolution between a host and a partially overlapping parasite can be divided into three distinct classes. First, when both species overlap only slightly, such that the coevolutionary hot spot is small, evolutionary dynamics are dominated by the selection experienced by host populations outside the geographical range of the parasite. Under these conditions, both species are monomorphic and the parasite species evolves to be perfectly adapted to its host, while the host is completely prevented from adapting to the parasite because of gene flow (Fig. 2, left hand regions). Secondly, as the size of the coevolutionary hot spot increases past a threshold value, stable allele frequency clines and spatial patterns of maladaptation and local maladaptation evolve in both species (Fig. 2, central regions). It is at this point, that the host begins to adapt to the parasite, thus causing a dramatic increase in parasite maladaptation. Thirdly, above a further threshold of hot spot size, synchronized allele frequency oscillations propagate throughout the geographical range of the host (Fig. 2, right hand regions). These synchronous oscillations preclude local maladaptation in either species.

Our results also demonstrate that the explicit location of the parasite species modulates the transitions between these different classes of dynamics. Parasites that overlap at the centre of the host's geographical range must do so to a greater extent for selection by parasites on hosts to destabilize the monomorphic equilibrium favoured by directional selection (compare Fig 2a, b, left hand lines). At the same time, however, centrally located parasites cause oscillations to evolve with significantly less overlap than those located at the host's periphery (compare Fig 2a, b, right hand lines). A consequence of these shifts in transition zones is that, all else being equal, a much smaller region of parameter space leads to polymorphic clines and local maladaptation when parasites occur at the centre of their host's geographical range.

These results have important implications for one of the central predictions of the geographic mosaic theory: trait mismatching and local maladaptation (Thompson, 1994; Thompson et al., 2002). Specifically, our results allow us to quantify the region of parameter space over which trait mismatching and local maladaptation are expected to occur. This region is depicted by the central area of Fig. 2. Inspection of this figure leads to two immediate conclusions. First, both species are more likely to exhibit trait mismatching and local maladaptation when rates of host gene flow are low. Secondly, trait mismatching and local maladaptation occurs more readily in both species as coevolutionary hot spots move toward the periphery of the host's geographical range (Compare Fig. 2a and b). Together, these results suggest that the primary predictions of the geographic mosaic theory will be manifested over the broadest range of hot spot sizes when rates of host gene flow are quite low.

Our results also suggest a novel mechanism that may generate observed empirical patterns of local maladaptation. In comparison with previous theoretical studies that have demonstrated local maladaptation because of asynchronous cycles (Gandon et al., 1996; Morand et al., 1996; Lively, 1999; Gandon, 2002), our results suggest that local maladaptation may also result from static allele frequency clines. Specifically, our results suggest that when coevolutionary hot and cold spots are coupled through the action of gene flow, permanent local maladaptation can arise only through stable clinal patterns. This mechanism may be particularly important given the strong tendency for coevolutionary cycles to synchronize in response to gene flow (Nuismer et al., 1999, 2000; Gandon, 2002).

Two unique predictions of the current model are amenable to comparison with observed empirical patterns. First, our model predicts that static clines in both host and parasite should be common when host and parasite overlap only partially. Secondly, our model predicts that patterns of local parasite adaptation should be constant across years. Both of these predictions receive at least some empirical support. The existence of static patterns of local adaptation receives support from studies of the parasite Microphallus and its snail host Potamopyrgus antipodarum. In an initial study of this interaction, Lively & Jokela (1996) found evidence for clinal variation in local adaptation. When this study was repeated one year later, Krist et al. (2000) found a very similar cline. Thus these two studies, that sampled from two different points in time, found similar clinal patterns of local adaptation. It seems unlikely that this result is compatible with the asynchronous cycles predicted by other models (Gandon et al., 1996). Similarly, clines of resistance and virulence have also been observed in the Drosophila parasitoid Asobara tabida (Kraaijeveld & Godfray, 1999).

In contrast, results from other studies seem to be more compatible with models that permit asynchronous oscillations. In a study of local adaptation in the anther smut fungus Microbotryum violaceum, for instance, Kaltz et al. (1999) found no evidence for organized clinal patterns of local adaptation. Similarly, Thrall et al. (2002) found no evidence for clinal patterns of local adaptation in the interaction between Linum marginale and its pathogen Melampsora lini. In a very different system, Little & Ebert (2001) found that patterns of local adaptation between Daphnia magna and Pasteuria ramosa shifted over time in two of three populations studied. These conflicting results seem to suggest that nuances of life history, population structure, or the genetic basis of coevolution may lead to patterns different from those predicted by the current model. Indeed, our model makes several genetic assumptions that may limit its applicability to a subset of host–parasite interactions. For instance, we assume that coevolution is governed by matching alleles at a single haploid locus. Consideration of additional loci (e.g. Seger, 1988; Peters & Lively, 1999), diploidy (Switkes & Moody, 2001), or interactions mediated by quantitative traits or gene-for-gene mechanisms (e.g. Doebeli, 1997; Gavrilets, 1997; Abrams, 2000; Sasaki, 2000; Agrawal & Lively, 2002) may lead to different quantitative conclusions. We anticipate, however, that the majority of our qualitative conclusions will be robust across a wide range of coevolutionary genetic systems.

A more crucial issue may be the extent to which our assumptions about the spatial arrangement of selection accurately reflect natural host–parasite interactions. Our model makes a very specific assumption that one of the alleles at the host interaction locus is consistently favoured by directional selection. This could clearly be the case if one allele had an associated cost or pleiotropic effect. While there is evidence that certain alleles or phenotypes responsible for host resistance may impose costs or have pleiotropic effects (e.g. Ruwende et al., 1995; Webster & Woolhouse, 1999), it is not clear that this is generally the case (Bergelson et al., 2001). If in some situations loci responsible for parasite resistance have no effect outside the geographical range of the parasite, our conclusions would no longer apply.

In a broader sense, our results may have consequences for the evolution of parasite geographical range. In particular, our results show that parasites will, in general, be best adapted to their host species when they occupy only a small fraction of the host's geographical range. This high level of adaptation can be achieved because host gene flow from populations outside of the parasite's geographical range effectively swamps a host's ability to adapt to the parasite. If, however, the geographical range of the parasite increases past some critical value, hosts overcome the swamping influence of gene flow and rapidly adapt to the parasite. Thus insofar as adaptation to a particular host limits the geographical range of a parasite, this mechanism may play a role in shaping the geographical ranges of parasites, relative to their hosts. This result builds on other recent results that have begun to suggest that gene flow may play a fundamental role in shaping the geographical ranges of species (Garcia-Ramos & Kirkpatrick, 1997; Holt & Gomulkiewicz, 1997; Kirkpatrick & Barton, 1997; Case & Taper, 2000; Nuismer & Kirkpatrick, 2003).

We have considered a simple model of coevolution between a host and a partially overlapping parasite. Despite its simplicity, our model has produced several novel and testable predictions for the expected spatial structure of allele frequency clines and patterns of local maladaptation. Although they are not rigorous tests of this model, the results of several empirical studies are compatible with these predictions. Alternative models that include more realistic genetics and explicit demography will be crucial to refining these predictions. Such models will also be necessary for the continued development of a predictive geographic mosaic theory of coevolution.

Acknowledgments

This work was supported by NIH F32 GM65620-01 to S.N., NSF DEB-0073911 to J.N.T and NSF DEB-9528602 to R.G. We thank Mark Kirkpatrick, Dave Hall, Martin Morgan, Daphne Fairbairn, and two anonymous reviewers for helpful suggestions.