## Abstract

Infections by multiple genotypes are common in nature and are known to select for higher levels of virulence for some parasites. When parasites produce public goods (PGs) within the host, such co-infections have been predicted to select for lower levels of virulence. However, this prediction is based on simplifying assumptions regarding epidemiological feedbacks on the multiplicity of infections (MOI). Here, we analyse the case of parasites producing a PG (for example, siderophore-producing bacteria) using a nested model that ties together within-host and epidemiological processes. We find that the prediction that co-infection should select for less virulent strains for PG-producing parasites is only valid if both parasite transmission and virulence are linear functions of parasite density. If there is a trade-off relationship such that virulence increases more rapidly than transmission, or if virulence also depends on the total amount of PGs produced, then more complex relationships between virulence and the MOI are predicted. Our results reveal that explicitly taking into account the distribution of parasite strains among hosts could help better understand the selective pressures faced by parasites at the population level.

## 1. Introduction

Host co-infections by several genotypes of the same parasite species are ubiquitous, and are found in several human infections [1]. Epidemiological models show that co-infections by different genotypes can select for more virulent strains in the host population provided that virulence confers a within-host competitive advantage to the parasite [2–5]. Such a correlation between the virulence of a strain and its within-host competitiveness has been demonstrated experimentally for parasites such as *Plasmodium chabaudi* in mice [6] or microsporidia in *Daphnia* [7].

The generality of the result that co-infections select for more virulent strains has been challenged by theoretical and experimental studies based on parasites producing public goods (PGs). The rationale is that, if high levels of parasite cooperation translate into high levels of virulence, then diverse infections should favour ‘cheating’ strains that do not cooperate, and therefore lower levels of virulence should evolve. Support for this idea has come from both empirical [8] and theoretical models of ‘collective action’ between viruses [9,10].

Many prokaryote species can synthesize siderophores, which are diffusible iron-chelating compounds that bind to iron with a high affinity and allow iron to be imported into the cell [11]. Because siderophores are costly to produce but free to use, bacterial strains that cheat and do not produce siderophores should have a fitness advantage when competing with siderophore-producing strains, and it should not pay to cooperate and produce siderophores. This scenario has been studied theoretically by West & Buckling [12]. Using a kin selection model in which siderophore production is positively linked to parasite virulence, they find that the optimal level of siderophore production decreases with the number of bacterial strains per host. In other words, they predict that low within-host parasite relatedness selects for lower levels of parasite cooperation (and therefore lower levels of virulence). This pattern was found experimentally by Harrison *et al.* [13] for the bacterium *Pseudomonas aeruginosa*.

Based on these studies, it has been suggested that co-infections generally select for less virulent strains when host exploitation depends on cooperation by parasites [10,14]. Here, we argue that the situation can be more complex than previously expected because of epidemiological processes. Indeed, in natural populations, the optimal level of virulence will be determined by the interplay between the selective pressures at both the host and population levels [3,15]. The two need not necessarily agree. First, natural populations are characterized by a heterogeneous distribution of parasite strains (some hosts are uninfected, others are singly or doubly infected, and so on), which has been shown to be shaped by parasite traits (e.g. a virulent strain that causes short infections should yield a distribution with few co-infections [16]) and to alter the selective pressures experienced by parasites (e.g. different strains can be selected for depending on the frequency of susceptible hosts in the population [17]). Second, evolutionary trade-offs at the epidemiological level have a crucial effect on virulence evolution. For instance, theory predicts that a saturating trade-off between transmission rate and virulence should lead to the existence of an intermediate level of virulence [3,5,18–20]. Such a trade-off has been documented empirically (for instance for HIV in humans [21] and for a protozoan parasite of monarch butterflies [22]). In all these cases, intermediate levels of virulence maximize parasite fitness.

To explore in more detail the evolutionary consequences of epidemiological feedbacks on the evolution of PG production by parasites, we use a nested model [23] (i.e. we analyse explicitly the interplay between within-host and epidemiological dynamics). This allows us to specify the relationship between parasite density and epidemiological variables (transmission rate and virulence), and to investigate the consequences of an uneven distribution of strains across hosts. Our within-host dynamics are inspired by the example of siderophore-producing bacteria, and are directly comparable to previous models (e.g. [12]), although our general approach and message are not restricted to this particular example. Between-host dynamics are based on a co-infection model that allows coexistence between up to two strains in a host [3]. By varying the multiplicity of infection (MOI) in the host population from 1 (only single infections) to 2 (only co-infections), we show that the result that the optimal virulence increases with the MOI is only valid if there is a linear trade-off relationship between virulence and transmission. When there is a saturating trade-off or when virulence depends on the total amount of PGs produced in the host, different patterns between virulence and MOI are predicted, which depend on the feedback between MOI and epidemiology.

## 2. Within-host dynamics

We consider a host infected by *n* parasite strains. The within-host density of strain *i* is denoted *N*_{i}. A strain is defined by its trait *x*_{i}, which corresponds to the units of PG it produces. The dynamics of the density *G* of PG is given by
2.1where *u* is the uptake rate of the PG. This equation assumes that all strains have the same uptake rate, and that the uptake of PG by strain *i* is proportional to its frequency *N*_{i} /*N*, where *N* = ∑_{i=1}^{n} *N*_{i} is the total within-host density of parasites.

The dynamics of strain *i* follows an ordinary differential equation with three terms reflecting, respectively, a logistic growth, a cost of PG production (which is proportional to the strategy of the strain, *x*_{i}) and a benefit of PG exploitation (which is a function of the uptake rate of PG by strain *i*). This yields
2.2where *K*/*n* is the maximum population size achievable by each of the *n* strains in the absence of PG, *c* is the cost of producing one unit of PG and *b* is the benefit of PG exploitation (figure 1*a*).

Our within-host model is directly comparable to earlier studies, in particular the siderophore production model of West & Buckling [12]. The main difference is that we work directly with within-host parasite population sizes, whereas West & Buckling [12] postulated a relationship between the within-host growth rate of parasites and PG production and exploitation. One advantage of our approach is that it leads to predictions in terms of parasite densities, which are in principle easier to assess than parasite growth rate. Note that, in order to allow for long-term parasite coexistence within the host, we have assumed that all strains are equally competitive to obtain baseline resources in the host (meaning that when there is no PG production, each strain can reach a maximum density of *K*/*n*). This assumption also allows us to distinguish between competition for host resources and competition for PG.

Following earlier studies [12,24], the amount of PG received by an individual parasite, *uG*/*N*, is scaled by a parameter *a* ≤ 1 to take into account the fact that the cost of producing the PG may increase faster than the benefit gained by the exploitation of the PG. As a consequence, when *a* < 1, there exists a PG production rate that maximizes parasite population size (figure 1*b*). In equation (2.2), we assume that the benefits and costs associated with PG uptake and production are additive, but our results also hold when they are multiplicative, as shown in the electronic supplementary material, S7.

We further assume that within-host dynamics are fast compared with between-host dynamics so that within-host densities can be assumed to be at equilibrium. The equilibrium density of PG is 2.3and the equilibrium density can be found by solving 2.4where 2.5is the average rate of PG production in the population, with . For two strains, there are two stable equilibria: one where both strains densities are zero and one where they are strictly positive. In the following, we restrict our attention to the range of parameter values such that within-host parasite densities are non-zero.

An important parameter in kin selection models of within-host evolution is the relatedness between parasite strains within a host. Clearly, relatedness will depend on the MOI, because it should be 1 in a host infected by a single parasite strain and less than 1 in a multiply infected host. Following Taylor & Frank [25], within-host relatedness can be defined more precisely as the regression of the average rate of PG production over the rate of PG production *x*_{m} of a focal strain, which yields for our model
2.6

Hence, if all strains are equally abundant, or (in our model) if all strains have the same trait value, we have a simple relationship between within-host relatedness and the MOI within a host, *r*_{m} = 1/*n*. In general, however, relatedness will also depend on how within-host dynamics shape parasite densities.

## 3. Epidemiological parameters

In general, the epidemiological parameters of the disease (i.e. virulence and transmission) will depend on the within-host densities of each strain [23] and on the within-host density of PG. Here, we explore two ways to link the within- and between-host models.

### (a) Model 1: transmission–virulence trade-off through within-host densities

Following the majority of previous nested models [23], we assume that the transmission rate of a strain (*β*_{i}) is linearly correlated with the within-host density of this strain, and that virulence (*α*), which is defined here as the increase in host mortality owing to the infection, is a power function of the total parasite within-host density. This gives
3.1and
3.2where *β*_{0} and *α*_{0} are proportionality constants, and *ξ* describes the shape of the transmission–virulence trade-off curve [18–20]. Note that since virulence depends on the total parasite density, it is not a linear function of the amount of PG produced (figure 1*b*).

If *ξ* = 1, transmission is a linear function of virulence and theory predicts that the parasite strain with the highest virulence will always take over the host population. If *ξ* > 1, transmission is a saturating function of virulence, which implies that there exists an evolutionarily stable strategy (ESS) for the parasite corresponding to an intermediate level of virulence [3,5]. In the electronic supplementary material, S4.3, we show that the predictions of our model are qualitatively similar whether transmission or virulence is assumed to be a nonlinear function of within-host parasite density.

### (b) Model 2: detrimental public goods

In model 1, we have assumed that virulence depends only on parasite density, but other relationships are possible. For instance, virulence could also be affected by the total density of PGs (e.g. if the PG is toxic, or if, as in the case of siderophores, it takes away resources from the host). We model this scenario with the following definition for virulence: 3.3

Thus, virulence is a function of the total parasite density () and of the total amount of siderophores produced (). We do not assume that virulence depends only on the amount of siderophores produced because a strain that does not produce siderophores would then be completely avirulent.

## 4. From epidemiological dynamics to parasite fitness

We consider a host population at equilibrium infected by a resident parasite strain with trait *x*_{w}. A new mutant strain with trait *x*_{m} will increase in frequency if its epidemiological fitness, given by the basic reproductive ratio *R*_{0} (*x*_{w}, *x*_{m}) [26], is greater than 1. Note that this fitness depends on the mutant trait and on the resident trait. If mutations have a small phenotypic effect on the trait, the direction of selection is given by the selection gradient ∂*R*_{0}/∂*x*_{m} evaluated at *x*_{m} = *x*_{w}. Following classical invasion analyses [27,28], potential evolutionary endpoints *x** can be found by solving
4.1

Our aim is to investigate the evolution of within-host cooperation (i.e. the value of *x*) and its consequences for epidemiological parameters (particularly virulence).

### (a) Epidemiological dynamics

We consider an epidemiological co-infection model introduced by van Baalen & Sabelis [3], with one additional twist. Following Alizon & van Baalen [29], we introduce a parameter *σ*_{I} that measures the susceptibility of infected hosts to the disease. When *σ*_{I} = 0, infected hosts are immune to co-infection, and we recover the classical susceptible–infected (SI) model with single infections only. In the limit where *σ*_{I} is very large, infected hosts immediately enter the co-infected stage, and the population is entirely composed of susceptible and co-infected hosts. For completeness, we also introduce a parameter *σ*_{S} that measures the susceptibility of uninfected hosts. The originality of this model is that it takes into account co-infections by the same strain, which is necessary to correctly assess the invasion fitness of a rare parasite strain (for further details, see [3,30]). We further assume that the order of infection does not matter (*D*_{mw} is equivalent to *D*_{wm}). This is in line with our assumption that recovery does not occur (i.e. infections are persistent). In other words, transient dynamics within a host (when a new strain arrives) can be neglected because they are short compared with the duration of infection. The life cycle is described in figure 2, and the full dynamics are given in the electronic supplementary material, S2.1.

The invasion fitness of a mutant strain [3] (electronic supplementary material, S2.2) is
4.2where and indicate the equilibrium densities of susceptible and infected hosts before the mutant appears. *λ*_{w} and *λ*_{m} are the force of infection of the resident and of the mutant strain, respectively (their full expressions are shown in the electronic supplementary material, S2.1), *α*_{wm} is the overall virulence in a host co-infected by the resident and the mutant strain, *β*_{mw} (respectively, *β*_{wm}) is the transmission rate of the mutant (respectively, resident) strain in a host co-infected by the resident and the mutant strain.

There are two interesting limit cases of equation (4.2). First, when *σ*_{I} = 0, we recover the classical expression of the basic reproductive ratio of a rare mutant parasite for an SI model with no co-infections [26] as
4.3

Second, when all hosts are infected in the population (), or when *σ*_{I} is very large, the evolutionary outcome will be determined solely by the ratio between parasite transmission and host mortality:
4.4In equation (4.4), parasite fitness at the between-host level is the product of the transmission rate of the focal mutant lineage and of host survival. This is similar to fitness expressions used in some previous kin-selection models [5,12], which we recover as a limit case of our model.

The comparison among equations (4.2)–(4.4) shows that, for a given trade-off between transmission and virulence, our model provides an additional level of epidemiological realism by explicitly taking into account the distribution of parasite clones among hosts. In our model, this distribution is characterized by the densities of uninfected, singly infected and doubly infected hosts, but generalizations to more complex distributions are possible. Finally, note that equation (4.2) does not depend on the host demography in the absence of infection, but only on the epidemiological processes of primary and secondary infections. In the remainder of the paper, we focus on a simplified model of host demography and assume that the total host population is constant (i.e. that every death is replaced by a new susceptible).

### (b) Feedback between epidemiology and the multiplicity of infections

Variations in the susceptibility of infected hosts to secondary infections, *σ*_{I}, will affect the average MOI in the population, which can be measured as
4.5where *I* and *D* are the total densities of singly and doubly infected hosts. In our model, the MOI takes values between 1 (when *D* = 0, i.e. all infected hosts are singly infected) and 2 (when *I* = 0, i.e. all infected hosts are co-infected). The quantity 1/MOI is the average within-host relatedness in the population and takes values between 1/2 and 1.

The precise relationship between MOI, *σ*_{I} and other epidemiological parameters will depend on the host demography in the absence of infection. For a population with a constant total size (*S* = 1 − *I* − *D*), we find at epidemiological equilibrium
4.6where is the equilibrium density of susceptible hosts.

Equation (4.6) shows that epidemiology affects the MOI in the population through the susceptibility to secondary infections *σ*_{I} and the density of susceptible hosts. In particular, we see that the MOI is an increasing function of *σ*_{I}, and a decreasing function of *S*. When the density of susceptible hosts is low, the MOI is expected to be close to 2. When *S* increases as a result of variations in demographic (*μ*) and epidemiological traits (*α* and *β*), so does the proportion of singly infected hosts. Therefore, ultimately, there is a feedback between investment into PG and the MOI in the population. This feedback will be affected by (i) the demography of the host population, (ii) the trade-off between virulence and transmission, and (iii) the rate of secondary infections. Taken together, equations (4.2) and (4.6) show that different mortality (varying *μ*) or infection (varying *σ*_{I}) regimes may have an impact on both the optimal PG investment and the average MOI in the population.

### (c) Evolutionary correlations between multiplicity of infection, public goods production and virulence

#### (i) Model 1: linear transmission–virulence trade-off

When the trade-off between transmission and virulence is linear (*ξ* = 1) and co-infections do not occur (*σ*_{I} = 0), selection favours strategies of host exploitation that maximize virulence [18,26]. In our model, this implies that PG production evolves until within-host parasite density is maximal (which occurs at *x* = 1 with our default parameter values; see figure 1*b*). When infected individuals become more susceptible to co-infection (*σ*_{I} increases), we find that the optimal level of PG production decreases (figure 3*a*). Because the average MOI is an increasing function of *σ*_{I}, this also results in a decreasing relationship between MOI and optimal PG production (figure 3*b*). The interpretation of this pattern is that singly infected hosts provide a refuge for cooperative strains, where they cannot be exploited by cheaters. When the proportion of singly infected hosts in the population decreases, this advantage weakens and cheaters are increasingly favoured. As a result, within-host parasite density decreases when the MOI increases, and so does evolutionarily stable (ES) parasite virulence (figure 3*c*). This corroborates the findings of previous models (e.g. [12]).

#### (ii) Model 1: saturating transmission–virulence trade-off

Assuming a saturating transmission–virulence trade-off (*ξ* = 2; figure 3*d*–*f*) leads to a radically different outcome. For high values of the susceptibility of infected individuals to secondary infections (*σ*_{I}), we find that PG production is never selected for. In this case, there is no correlation between virulence and MOI. By contrast, for low values of *σ*_{I}, evolutionary bistability is observed. If the population starts from a low investment into PG production, PG production is always counter-selected and the MOI has no effect on virulence evolution. If the initial PG investment is above a threshold, selection favours high investment into PG production, and this optimal level of PG production decreases as infected hosts become more susceptible to secondary infections. In the limit where only single infections are possible (*σ*_{I} = 0), we show in the electronic supplementary material, S4 that, if 0 < *a* < 1, the threshold *x*_{c} above which PG production is favoured is *x*_{c} = (*c*/(*ba*))^{1/(a−1)}, which yields *x*_{c} = 1 with our parameter values. Hence, two evolutionary outcomes are selected for depending on the initial investment in cooperation: either pure defection (no PG production) or an intermediate investment into PG production.

Understanding what happens in the limit where only single infections occur sheds light on the processes underpinning this result. For a saturating trade-off, selection favours intermediate levels of virulence. Hence, strains maximizing within-host parasite density (*x* = 1) are counter-selected, and this strategy actually becomes an evolutionary repellor. Because virulence is a bell-shaped function of PG production, an intermediate level of virulence can be achieved through either a low or a large investment into PG production, depending on which side of the repellor the population starts off (figure 4). Now, as the susceptibility to secondary infections increases, so does the advantage to cheaters. This has two effects. First, a zero investment into PG production becomes attractive for a wider range of initial conditions (the repellor moves to the right in figure 4*b*). Second, the non-zero ESS decreases because it pays less to be a cooperator, which increases virulence (the ESS moves to the left in figure 4). At a critical value of *σ*_{I}, the two points collide and 0 becomes the sole evolutionary endpoint.

Because the MOI and the optimal PG production both depend on the susceptibility parameter *σ*_{I}, the correlations between MOI and virulence are more complex under this scenario (figure 3*f*). Our model predicts either no effect of MOI on virulence, or that virulence increases with average MOI despite selection for lower PG production rates. In particular, populations can display identical levels of MOI, but very different levels of virulence, depending on the value of *σ*_{I}. Hence, different epidemiological dynamics can lead to similar levels of MOI at evolutionary equilibrium, but opposite predictions for virulence evolution.

In the limit where *σ*_{I} → ∞, our model converges to the limiting case where all hosts are co-infected. This is the scenario investigated by previous models (e.g. [12]). These models used a linear trade-off and thus predicted an intermediate optimal PG production, but we further show that the evolutionary outcome actually depends on the shape of the epidemiological trade-off. In the electronic supplementary material, S5.1, we show that if *ξ* > 1 + *μ* /*α*_{0} (i.e. if the trade-off is sufficiently concave), siderophore production cannot evolve from a low initial investment into siderophore production. If *ξ* is below this threshold value (in particular if the trade-off is linear), selection selects for a low, but non-zero investment in PG production, as predicted by West & Buckling [12]. In the electronic supplementary material, S6, we further show that this result holds qualitatively for the whole range of susceptibility values. Note that the predictions of a linear trade-off are not robust for even slightly nonlinear trade-offs (e.g. *ξ* ≈ 1.1).

#### (iii) Model 2: detrimental public goods

We now turn to model 2 and assume that parasite virulence depends both on within-host parasite density and on the total amount of PG produced (equation (3.3)). We find a non-monotonous relationship between the optimal level of PG production and the susceptibility to secondary infections (figure 3*g*). As a result, selection leads to both PG production and virulence being maximal at intermediate MOI (figure 3*h*,*i*). The reason why this happens is that the direct contribution of PG to virulence yields an additional cost in terms of parasite survival (compared with the scenario where the transmission–virulence trade-off is linear, virulence increases faster). In the electronic supplementary material, S4.2, we show that this additional cost becomes stronger when the MOI decreases, whereas the benefits of cooperating go in the opposite direction. This explains why the net benefit-to-cost ratio of cooperators peaks at an intermediate MOI.

#### (iv) Effect of background host mortality

Equation (4.6) shows that the MOI is also a function of other epidemiological traits such as host mortality *μ*. Everything else being constant, the MOI is a decreasing function of *μ*. In the electronic supplementary material, figure S4, we show the resulting correlations between MOI, PG production and virulence when varying host background mortality and keeping the susceptibility to secondary infections *σ*_{I} constant. Comparing figure 3 with electronic supplementary material, figure S4, reveals that different evolutionary correlations may be expected between MOI and virulence in the population depending on which epidemiological parameter varies across populations or treatments. Although the prediction that virulence decreases with the MOI when the trade-off is linear is qualitatively preserved, variations in host mortality or in susceptibility do not in general lead to the same patterns for other trade-off functions. Note that the complexity of this result is not surprising since the value of *μ* affects at the same time the optimal PG production rate within a host (by determining the duration of the infection) and the MOI (which itself affects the optimal PG production through epidemiological feedbacks). This suggests that, to fully understand the relationship between virulence and MOI at the population level, a detailed knowledge of the mechanisms underpinning the change in MOI is required.

## 5. Discussion

The current paradigm is that, for parasites producing PGs in their hosts, such as phages [8] or siderophore-producing bacteria [12,31], multiple infections should select for lower levels of virulence [10,14]. Indeed, if virulence is associated with the amount of PGs parasites produce, co-infections can favour less virulent (‘cheater’) strains because they allow them to exploit more virulent (‘cooperative’) strains. In other words, when within-host relatedness decreases, parasite cooperation should decrease, and so should virulence [9,12,13,32]. Here, we show that these predictions need not always hold because epidemiological feedbacks and host demography can alter the selective pressures on virulence through changes in the distribution of parasite strains among hosts.

Using a mechanistic within-host model of PG production embedded in an epidemiological framework allowing for co-infections, we explore the evolutionary correlations between the MOI and parasite virulence under different assumptions regarding the relationship between virulence, transmission, PG production and parasite within-host densities.

With a linear trade-off between transmission and virulence (as assumed by West & Buckling [12]), selection favours ever-increasing levels of virulence, and therefore population-level feedbacks do not constrain within-host investment in PG production. Hence, when only single infections are permitted, the ESS is the level of PG production that maximizes within-host parasite density. By contrast, when the fraction of multiple infections among infected hosts increases, the optimal investment in PG production will be lower because defectors have an advantage in co-infected hosts.

The evolutionary outcomes are different when we assume a saturating transmission–virulence trade-off, which is at the core of many models in evolutionary epidemiology [20]. In particular, the optimal level of virulence can be higher when the MOI is maximal than when it is minimal, which contrasts with earlier predictions. This is due to the fact that, when virulence increases faster than transmission, population-level feedbacks favour intermediate levels of virulence because when virulence is too high, parasites do not transmit sufficiently before the host dies [18–20]. Here, such intermediate levels of virulence are achieved when within-host parasite density is not maximal, which occurs at both low and high (but not intermediate) levels of PG production. When the MOI in the population is sufficiently high, only a zero investment in PG production is evolutionarily stable, and there is no correlation between MOI and virulence. By contrast, when the MOI in the population is low to moderate, assuming a saturating trade-off introduces a bistability, and non-zero levels of PG production can evolve if the initial level of PG production is above a threshold. Otherwise, PG production is counter-selected. This results in a negative relationship between MOI and PG production, but a positive relationship between MOI and virulence.

In general, the exact shape of the trade-off between transmission and virulence will affect evolutionary predictions. For instance, assuming that virulence depends both on parasite and on PG within-host densities, we find that PG production and virulence can peak at intermediate levels of the MOI, resulting in a bell-shaped relationship between MOI, PG production and virulence. A possibility for a sequel study would be to explore trade-off shapes more extensively by resorting to a critical function analysis [33].

Our analysis shows that the relationship between the MOI, PG production and parasite virulence may be less general than expected from previous studies that showed that, for PG-producing parasites, multiple infections should always select for lower PG production and lower levels of virulence. There are two main reasons for this. First, some studies (e.g. [12]) assumed a linear transmission–virulence trade-off, which is unrealistic for many infectious diseases and leads to more virulent strains always having an advantage at the between-host level. Hence, when the trade-off is linear, virulence evolution is only constrained by within-host parasite interactions, while parasite between-host fitness is known to be shaped by population-level processes that can alter and, as we show, even reverse within-host selective pressures [26]. Second, the MOI in previous studies (reviewed by Brown *et al.* [10] and Buckling & Brockhurst [14]) was usually defined as a fixed parameter *n* that can be varied freely (‘open models’ [34]), whereas what empiricists generally have access to is a population average that will be affected by epidemiological and demographic processes (e.g. the MOI generally depends on transmission, virulence, and other host and disease life-history traits).

### (a) Kin selection

In our study, the MOI is measured as the population average of the within-host relatedness between parasites, and is therefore affected by epidemiological and demographic feedbacks. It is therefore not surprising that our study has strong links with demographic kin selection models of the evolution of social traits. Indeed, an alternative way to view our epidemiological setting would be to consider the population of hosts as an infinite island model [35,36], each host being an island that can be infected by 0, 1, … , *n* parasite strains. This is analogous to Alizon & Taylor's [37] model of fecundity helping, with the important differences that in our model parasite clones play the role of individuals, and transition rates are affected by within-host densities, not only by the number of clones sharing a host. Our model is actually a limiting case of this process where we set *n* = 2 for simplicity. Indeed, this allows us to derive a simple analytical expression for the *R*_{0} of a mutant parasite. Moreover, previous studies have also shown that models with at most two strains per host already capture much of the effect of multiple infections on virulence evolution [3,29,38].

This analogy allows us to clarify the key difference between our model and previous studies of PG production in parasites, which typically consider that mutations have no effect on the distribution of parasite strains across individuals. This amounts to neglecting selection through the effects of the trait on the demographic state of the population, which can be an important factor influencing the evolution of social traits such as cooperation and virulence [39]. By contrast, our model allows us to consider how such demographic or epidemiological processes affect the distribution of parasite strains across hosts, and how the resulting variation in MOI alters selection on PG production and virulence.

### (b) Empirical implications

Harrison *et al.* [13] have shown experimentally that, when co-infecting a host with two strains with different levels of virulence, the less virulent strain was favoured. However, their experimental set-up follows closely the model assumptions of West & Buckling [12] and focuses on isolated hosts. Our analysis suggests that to fully understand how multiple infections affect the evolution of virulence, population-level processes (such as transmission between hosts) need to be considered in future studies. This has been done in other empirical studies that take into account demographic and epidemiological feedbacks (e.g. [40,41]).

A further empirical implication of our results is that the relationship between PG production and virulence is not necessarily linear (figure 4). Hence, comparing two strains with different levels of virulence may not be sufficient to analyse the selective pressures on virulence in PG-producing parasites and to determine the evolutionary outcome. To gain more insight, future experiments should first test the link between PG production rate and host virulence. This could be achieved by generating a collection of bacterial strains that vary in the siderophore production rate (as in [42]) and assessing parasite density and virulence for each strain. This would allow empiricists to test the bell-shaped within-host trade-off between parasite density and PG production that is typically assumed in theoretical models, and also shed light on the range of biologically feasible levels of PG production. Our model reveals that, owing to epidemiological feedbacks, selection could favour overproduction of PG beyond the value that maximizes within-host parasite density. An important question is whether such levels of PG production are biologically feasible. Arguably, this overshooting may be too costly at the within-host level and prevented by selection for regulatory mechanisms. Such mechanisms have been studied in *P. aeruginosa* and data show that regulation can be facultative, or lost in some mutants [43,44]. This calls for further empirical and theoretical investigations of the joint evolution of siderophore production and regulatory mechanisms.

Testing the effect of multiple infections on virulence evolution is notoriously difficult. One of the main challenges is that it is technically complicated to prevent multiple infections from occurring. A possibility would be to use the setting developed by Ebert & Mangin [40], which simulates a large host baseline mortality by removing a fraction of hosts and replacing it with the same amount of susceptible hosts. This marginalizes multiple infections but it can also generate complicated short-term evolutionary dynamics [45]. Here, we show that in addition to this issue, the way the MOI varies (e.g. host susceptibility to co-infections versus host background mortality rate) is likely to affect the evolutionary outcome of the model.

### (c) Perspectives

Few co-infection models consider within-host dynamics explicitly (but see [46]), which means that theoretical predictions often cannot be interpreted in terms of within-host cooperative or competitive interactions. An exception is provided by Alizon & van Baalen [29], who analysed a model of within-host competition between parasites mediated by the immune system. Here, we focus on PG-producing parasites, but there exists a huge variety of within-host interactions between parasites, which deserve further investigation [47].

Throughout most of our analysis, we limited the number of strains per host to two. Classical open models do not have this limitation, but they also ignore epidemiological feedbacks. Our simplification allows for analytical tractability, while capturing the important features of the co-infection process [3]. Moreover, because we assume rare mutations, evolution proceeds through the replacement of one strain by another and the parasite population is at most dimorphic. However, the interplay of mutation and infection processes could generate a higher between-host variance in the distribution of parasite strains. To take into account the fact that hosts can harbour more than two parasite genotypes, alternative evolutionary models would be required (e.g. models of co-infection by different species [38] or population-genetic models [45]). This would provide some interesting insights into the robustness of our result to a higher genetic variance and higher mutation rates within the parasite population.

Another simplifying assumption of our study is that the host population is well mixed. Yet the spatial or social structure of the host population has been shown to have a strong impact on disease evolution [5,48–50]. In a spatially structured host population, virulence evolution is known to be determined by a balance between epidemiological structuring and genetic structuring [50]. It will be interesting to investigate in more detail how the interplay between the within-host and between-host structures of the parasite population affects the evolution of virulence in PG-producing parasites.

## Acknowledgements

We thank S. Gandon for discussions, and T. Day, Y. Michalakis and three reviewers for helpful comments. S.A. is funded by the CNRS and the IRD. S.L. was supported by a Marie Curie Intra-European Fellowship FP7-220080 (Mobvir) from the Seventh Framework Programme of the European Commission, by CNRS and by an ANR grant RPDOC 2010 (SpatEvolEpid).

- Received March 9, 2011.
- Accepted April 15, 2011.

- This journal is © 2011 The Royal Society