## Abstract

We study the epidemiology of a viral disease with dose-dependent replication and transmission by nesting a differential-equation model of the within-host viral dynamics inside a between-host epidemiological model. We use two complementary approaches for nesting the models: an agent-based (AB) simulation and a mean-field approximation called the growth-matrix (GM) model. We find that although infection rates and predicted case loads are somewhat different between the AB and GM models, several epidemiological parameters, e.g. mean immunity in the population and mean dose received, behave similarly across the methods. Further, through a comparison of our dose-dependent replication model against two control models that uncouple dose-dependent replication from transmission, we find that host immunity in a population after an epidemic is qualitatively different than when transmission depends on time-varying viral abundances within hosts. These results show that within-host dynamics and viral dose should not be neglected in epidemiological models, and that the simpler GM approach to model nesting provides a reasonable tradeoff between model complexity and accuracy of results.

## 1. Introduction

Standard epidemiological models of viral infections assume that each host is infected with the same amount of virus and experiences the same time course of infection. Yet, in natural epidemics, both disease time courses and doses of transmitted virus vary among hosts. Some authors have attempted to integrate viral dose into epidemiological models, finding that it can determine invasion success and persistence of viruses (Dushoff 1996; Regoes *et al*. 2002; Li & Wang 2009; Yu *et al*. 2009). For example, Dushoff (1996) studied a modified susceptible–infected–recovered (SIR) model in which there are two classes of infected individuals: one infected with a high viral dose and more severe infection, and the other with a low viral dose and more mild infection. This modified SIR model predicts that a virus can persist in a population by means of transmission through the hosts infected with a large dose even if the virus would not persist in a regular SIR model. Furthermore, whether or not an infection can spread in a population when viral dose is considered depends on a threshold viral abundance independent of *R*_{0}, the average number of new infections caused by a single infected individual at the start of an epidemic (Regoes *et al*. 2002; Li & Wang 2009; Yu *et al*. 2009).

These results underscore that viral dose can determine epidemic dynamics. Yet, they have limited implications on understanding its impact on disease emergence, for two reasons. For one, dose is treated as a constant, which artificially amplifies its effects. Although constant dose size is an important extreme case to be considered for comparison, it is unlikely to occur commonly in nature. Secondly, in the epidemiological models of dose effects to date, dose-dependent dynamics within hosts are not considered. Therefore, these models cannot provide a mechanistic understanding of how viral dose impacts host immunity at the within-host and the population scale. To understand the role of dose- and time-varying viral dynamics on disease emergence risk, methods for analysing the interdependence of dose, within-host dynamics and transmission need to be developed.

The effect of viral dose on within-host dynamics *in vivo* is not completely understood, but some general principles have emerged. Viral dose impacts the viral growth curve as well as the strength of the host-immune response (Hughes *et al*. 2002; Quan *et al*. 2009). In viruses that cause acute infections, such as foot-and-mouth disease and influenza viruses, there is an inverse relationship between the size of viral dose and the infectious period, a critical determinant of transmission probability (Larson *et al*. 1976; Ottolini *et al*. 1996, 2005; Hughes *et al*. 2002; Alexandersen *et al*. 2003; Quan *et al*. 2004, 2009; Howey *et al*. 2008). Further, although the particular effects of viral dose on adaptive immunity levels are not completely understood, the data indicate that dose does in fact influence the dynamics of immunity, and that the relationship is complex and nonlinear. For example, the number of antibodies developed after clearing an infection increases with a very large viral dose (Ottolini *et al*. 1996; Hughes *et al*. 2002; Wherry *et al*. 2002), but a dose of moderate size tends to produce an amount of antibody comparable to that produced from a dose of small size. Thus, viral dose can influence both a central component of viral fitness, viral transmission and a major source of selection, host immunity.

Here, we develop two nested models (Mideo *et al*. 2008) for examining the effects of dose-dependent replication and transmission on epidemic behaviour. A single within-host model of viral and immunity dynamics is nested within two separate between-host models: an agent-based (AB) stochastic simulation model and a deterministic compartmental model that uses a next-generation matrix (GM) (Diekmann *et al*. 1990). These two between-host models differ in that the AB model is structured in the sense that each individual can only contact one other individual at each time step, whereas the GM model is mass-action, assuming an unhindered supply of connections between susceptible and infected individuals. A second difference is that the AB model allows for transmission during all possible states of the within-host dynamics, while the GM model coarse-grains within-host dynamics into discrete states that in turn determine transmission probability. To investigate whether the GM model is an accurate approximation to the more comprehensive AB model, we run simulations of epidemics under identical conditions for both models and compare model outputs. We show that the GM model is a reasonable approximation to the AB model. Further, through the use of two control models that remove dose-dependence on replication and transmission, we find that population immunity after an epidemic is qualitatively different when dose-dependent replication and transmission is considered.

## 2. Material and methods

### (a) Within-host model

We consider a virus that causes an acute infection. The virus replicates with a rate *r* to an intrinsic carrying capacity *K*_{V} and is removed by innate immunity *I*, non-specific memory cells *N* and specific memory cells *P*. Thus, we model three components of immunity based on the observed biology concerning acute viral infections such as influenza. Specifically, innate immune defences are the first to respond to viral infection, whereas non-specific memory cells begin being effective on the time scale of several days post infection (Hikono *et al*. 2006). With these principles in mind, we assume that *I* grows at a baseline rate of *b*_{I}(1 − *I*/*K*_{I}, where *b*_{I} is a constant rate of growth and *K*_{I} is the carrying capacity for *I*, and increases in proportion to viral abundance with rate constant *a*_{I}. On the other hand, production of *N* is delayed and grows in proportion to viral abundance with rate constant *a*_{N}. We model this delay using the Heaviside step function, * Θ*(

*t*).

*(*

*Θ**t*) = 0 for

*t*< 0 and

*(*

*Θ**t*) = 1 for

*t*≥ 0.

Specific memory cells arise from naive memory cells with a constant lag of *τ*_{P} and grow in proportion to the product of viral abundance and pre-existing *P*-cell abundance with rate constant *a*_{P}, and in proportion to naive memory cells with rate constant *c. N* decays with rate constant *d*_{N}, while *P* reaches an asymptote once *N* and *V* are cleared completely. We do not implement waning of specific memory cells since we do not consider reinfections here. Since the viral load does not drop until T cells are produced to fight the virus (Hikono *et al*. 2006), we assume that the innate immune response is less effective at lowering the viral titre than the T cell response, i.e. *k*_{I} < *k*_{N} < *k*_{P}, where *k*_{I}, *k*_{N} and *k*_{P} are the rate constants of removal of virus for *I*, *N* and *P*, respectively.

The equations governing the within-host dynamics are:

We chose the values of our parameters such that the relative values and the behaviour of our model conformed to the observed viral and immune dynamics. For example, the lag time for naive and specific immune cells has been set to 2.5 and 3 days, respectively, to fall in line with the observed progression of an immune response to an acute viral infection (Vermaelen *et al*. 2001; Legge & Braciale 2003; Lawrence & Braciale 2004; Hikono *et al*. 2006). We limit viral load to 10^{11} ml^{−1}, which agrees with observed maximum viral loads (Baccam *et al*. 2006; Stephenson *et al*. 2009). We set *k*_{P} to be considerably larger than either *k*_{I} or *k*_{N} to portray the superior efficacy of specific immune cells in clearing a viral infection (Cerwenka *et al*. 1999; Hikono *et al*. 2006). Note that effects of *P* on viral load are negligible in primary infections (*P*(0) = 0) because of the lag in *P* production, which captures the observed data (i.e. enhanced clearance from specific immunity in secondary infections but relatively little effect in primary infections). In primary infections, it is the levels of *P* that are produced from the viral dynamics that are the important focal point. The specific parameter settings for each simulation run are given in the captions of the corresponding figures.

### (b) Between-host models

We modelled transmission between hosts using two approaches, both of which incorporate our within-host model into an SIR-type compartmental model of the host population. Our first transmission model, the growth-matrix (GM) model, is completely deterministic. It condenses the continuous-time within-host viral dynamics into a distribution of transmission events, thus approximating the effects of viral infection dynamics on transmission during an epidemic. Our second model is an AB simulation that tracks the state of infected hosts throughout their infectious period and modifies transmission according to viral load at contact. The contact structure of hosts in the AB model is such that at each time step of fixed size, *Δ**t*, every host is randomly matched pair-wise with another host. Thus, transmission in the AB model is stochastic because of contact among particular individuals, and the GM model is a deterministic approximation of the AB model.

In order to effectively compare predictions of the two models, we make the following common assumptions. First, both models assume mass-action dynamics, a constant population size (*N* = 1000) and no reinfection of recovered hosts. Second, we assume that bottleneck size at transmission is constant (10^{−4}) such that the dose transmitted at any time *t*, is *V*(*t*) × 10^{−4}. We also set the probability of transmission to 1 when the viral load of the infectious host is greater than a minimum threshold of infection (1 virion), so that an infection does not successfully grow in a host unless the dose is greater than or equal to one viral particle. Since transmission can occur at each epidemic time step in the AB model, we set the transmission rate, *α* (see below), equal to 1 in the GM model.

### (c) Model 1: growth-matrix

To analyse the spread of virus through a host population that is heterogeneous with respect to within-host infection dynamics, we used a next-generation matrix (Diekmann *et al*. 1990), which operates in discrete time. The principle behind the next-generation matrix (*G*) is to divide the population into groups of individuals, classifying them by the dose received and, in the case of heterogeneous initial immunity, *P*(0) value. Given two classes of infected individuals, *i* and *j*, row *i* and column *j* of *G* represents the number of new infections of type *i* caused by one infection of type *j* in a single time step. We describe the potential doses and, in the case of heterogeneous initial immunity, *P*(0) values for infection of hosts by the vectors, *ϕ*(*t*) and *ϕ*_{R}(*t*), in which each entry describes the number of infected individuals who have received a certain dose with a given initial immunity, *P*(0).

### (d) Uniform initial immunity

All individuals start out susceptible with *P*(0) equal to some constant value across the population. We introduce one individual infected by a specific viral dose. This host has the same *P*(0) as all other hosts in the population. *ϕ* is then a vector in which the number of entries equals the number of doses we consider, where *ϕ*(0) has all entries equal to 0 except for the entry corresponding to the dose in the initial host which has a value 1. Our growth matrix, *G*, is constructed by the following formula:

Here, *β* is the rate of transmission. Recall that we have set *β* = 1 in order to agree with transmission in the AB model. *T*_{j} is the time at which 10^{−4} *V*(*t*) drops below the viral dose size corresponding to *i*, given the within-host dynamics for a viral dose size corresponding to *j*. If 10^{−4} *V*(*t*) never reaches the viral dose size corresponding to *i*, then we set *G*_{ij} equal to 0. *Δ**t* is the time step between updating the infection status of hosts and *n*_{j} is an integer count of the number of time steps at which 10^{−4} *V*(*t*) is equal to the viral dose size corresponding to *i*.

In a similar spirit, we construct a recovery matrix, *R*, which describes the movement of infected individuals into the recovered class. From the within-host dynamics, we can determine a time, *T*′_{J}, at which *V*(*t*) drops below 1 for a given viral dose *J*. We then have:

Here, *k* is the number of viral doses we consider.

Now, with a set population size, *S*(0), we can add *G**ϕ* to *ϕ*, subtract *R**ϕ* from *ϕ*, and add *R**ϕ* to *ϕ*_{R} repeatedly in that order until the sum of the entries in *ϕ* and *ϕ*_{R} is either just below or equal to *S*(0). We then apply *G* once more, but scale the new inflow of infected individuals so that the sum of the entries in *ϕ* and *ϕ*_{R} equals *S*(0). In this way, the number of new cases will be *G**ϕ* at each time step, and the last number of new cases will be a scaled value of *G**ϕ*. At this point, we can subtract *R**ϕ* from *ϕ* until the sum of the entries in *ϕ*_{R} reaches *S*(0) and the magnitude of *ϕ* drops below 1. The between-host dynamics are then given by the following equations:

Since we assume each individual that receives the same dose has the same course of infection, it is possible to calculate the mean immunity of newly recovered individuals simply by taking an average of the immune values for *R**ϕ* at each time step. Likewise, we can calculate the mean dose at each time step by taking an average of the dose received in *G**ϕ* at each time step. *R*_{0} is calculated as the ratio of the magnitudes of *G**ϕ* and *ϕ*. This calculation of *R*_{0} does not give the total number of new cases caused by an infected individual during the epidemic, but instead it produces a mean number of new cases in one time step caused by an infected individual.

### (e) Heterogeneous initial immunity

To consider heterogeneous *P*(0) in a population, we consider classes of infected individuals based on both the viral dose received and the value of *P*(0). Thus, for *m* levels of initial immunity and *k* viral doses, *ϕ* will have *mk* entries. To create a growth matrix for this situation, we assume that the number of new infections of a given dose caused by one individual can be divided evenly among the infected classes with any *P*(0) and that given dose. In other words,

Once again, we create a recovery matrix that is a diagonal matrix, such that the *i*th entry is the rate of recovery for the *i*th class.

In order to analysze infection spread in this case for a given initial dose, we run the same method described above for an initial host with 10 different values of *P*(0) equally spaced in [0, 2] and the *V*(0) considered. Epidemiological parameters are calculated in the same fashion as described above.

### (f) Agent-based model

We use an AB simulation model to examine whether condensing the within-host dynamics via the next-generation matrix approach is similar to epidemiological dynamics determined by the entire viral growth curve and random host contact. Each simulation begins with a single infected host. Ten replicate simulations are run for each set of initial viral dose, *V*(0), and specific immunity level, *P*(0). Contact between hosts occurs in discrete time. At each time step, half of the host population is randomly matched with the other half of the population (i.e. random contact of the entire population). The hosts are organized pair-wise such that each host has the opportunity to transmit the virus to exactly one other host. Transmission occurs only from an infected host to a susceptible host. For the case of heterogeneous initial immunity in the population, we assign a random *P*(0) value in [0, 2] to each host upon infection.

The number of infections caused during the infectious period of each host, *R*_{0}, is recorded at recovery. Mean *R*_{0} of recovering hosts is calculated at each time step. Likewise, mean dose in newly infected hosts, total new cases and mean peak-specific immunity levels in recovering hosts are calculated at each epidemic time step.

## 3. Results

### (a) Within-host model

The within-host model describes dose-dependent replication and immune response. Figure 1 shows the combined effects of viral dose and initial specific immunity level, *P*(0), on within-host behaviour. The maximum viral load and level of newly acquired specific immunity increases with viral dose for a given level of *P*(0) (figure 1*a*,*b*). On the other hand, both infectious period and time of maximum viral load decrease with viral dose (figure 1*c*,*d*).

A noteworthy aspect of our within-host system is the nonlinear relationship between *P*(0) and the strength of dose-dependent effects on replication dynamics. For example, both peak viral load and infectious period depend more strongly on viral dose at the extreme values of *P*(0), i.e. *P*(0) = 0 and 2, whereas only very small dose-dependent effects occur at moderate levels of *P*(0), i.e. *P*(0) = 0.5 and 1 (figure 1*a*,*d*). Similarly, *P*(0) has a different functional relationship to viral dose at moderate levels of *P*(0) relative to the extreme values, albeit generally producing the shortest time to peak viral load and highest level of peak immunity at moderate level of *P*(0) (figure 1*b*,*c*).

### (b) Population dynamics

In order to determine effects of within-host dynamics on epidemic behaviour in a population comprised entirely of susceptible hosts, we ran population level simulations using both the GM and AB models of transmission for a host population with no specific immunity (i.e. *P*(0) = 0 for all hosts). Owing to inherent differences in the approaches, there are some discrepancies between predictions from the two models. For example, the GM model predicts faster, more severe epidemics relative to the AB model, because of a methodological difference in the implementation of host contact. In the GM model, the infectious class has a constant supply of connections to susceptibles, whereas the AB model uses an explicit contact structure, which limits the rate of transmission since each host is connected pair-wise to only a single other host at each time step.

Despite this methodological difference in modelling contact, the models produce qualitatively similar and typical epidemic behaviour: initial exponential growth of new cases followed by a decline to extinction as the susceptible pool is depleted (figure 2*a*,*b*). Neither model shows a dependence of epidemic dynamics on dose in the first host, although the models do predict qualitatively similar effects of dose-dependent replication and transmission on mean immunity and dose trajectories. Mean immunity of recovered individuals starts out high but drops as the epidemic progresses (figure 2*c*,*d*). This decline corresponds to a decline in viral dose, which emerges as susceptible hosts become rare. Mean dose increases rapidly at the outset but then remains roughly constant throughout most of the epidemic, also until susceptible hosts become rare (figure 2*e*,*f*).

By contrast, the models predict slightly different trends in *R*_{0}, which is partly because of the different methods used for estimating *R*_{0} (figure 2*g*,*h*). While the AB model produces *R*_{0} based on the entire infection course of individual hosts that recover at a particular time step, the GM model produces an *R*_{0} based on the new cases produced in one time step. For this reason, the *R*_{0} produced by the GM model is lower. Further, a gradual decline in *R*_{0} occurs with the AB model but not in the GM model, a disagreement that can be attributed to the fact that the AB model captures the loss of frequent contact between susceptible and infected individuals as susceptibles become more rare, while the GM model does not.

### (c) Heterogeneous P(0)

In reality, the distribution of initial immunity in a population is almost always heterogeneous. For this reason, we ran simulations under heterogeneous population *P*(0). For the GM model, we introduced a new dimension to the classes of infected individuals in the form of a discrete number of *P*(0) values. For the AB model, *P*(0) for each host in the population was drawn randomly from a uniform distribution on the interval [0, 2]. For the GM model, the runs consisted of varying the initial host such that the dose remains constant, but *P*(0) was set to each of 10 levels considered. We assumed that *P*(0) in the population was uniformly distributed in the range of 10 values mentioned previously. For the AB model, the runs consisted of a constant viral dose in the initial host and a random initial immunity in the introductory host.

The models show qualitative agreement in their mean behaviour. Specifically, all replicate runs for both models reach a common equilibrium in mean immunity, mean dose and *R*_{0} (figure 3). This emphasizes that epidemiological behaviour depends on the within-host dynamics in currently infected individuals (rather than those of the introductory host), which is a common feature of our models. We also see that the mean immunity of recovered individuals gradually declines in both models (figure 3*c*,*d*), a result we noted also with uniform initial immunity in the population.

Both the number of new cases and the mean dose are notably higher in the GM model, and this can be attributed to the method used for modelling heterogeneous *P*(0) with the GM model. We assumed that upon infection, an individual will generate an equal number of new infections with a given dose size for each level of initial *P*(0) in the new host. Thus, completely susceptible individuals are always infected and are able to drive up the number of new cases and also the mean dose. Since new infections in the AB model have a stochastic element in the initial *P*(0) value, this effect from completely susceptible individuals is not present in the way it is for the GM model.

### (d) Effects of viral dose on epidemic behaviour

Next, we examined the effects of viral dose and initial immunity on epidemic behaviour and the profile of host population immunity via the following characteristics: (i) epidemic severity: the peak number of new cases; (ii) epidemic rate: the time to the peak number of cases (figure 4); (iii) mean population immunity: mean immunity in the host population after the epidemic; and (iv) variability in immunity: the coefficient of variation in this mean (figure 5). We note that the dose in the initial host does not have much effect in this analysis, except under heterogeneous *P*(0) in the AB model and *P*(0) = 2 in both the models. For heterogeneous *P*(0) in the AB model, we see that for a dose of 1 virion in the initial host, it is very likely that the epidemic will not spread. Similarly, for *P*(0) = 2 in both models, a large viral dose in the initial host is necessary for epidemic transmission.

The effects of initial immunity on epidemic severity and rate reveal some disagreement between the AB and GM models. First, in the AB model, epidemic severity is sensitive to moderate increases in *P*(0), whereas the GM model predicts similar case rates across low to moderate *P*(0) (figure 4*a*,*b*). Furthermore, at high *P*(0), a positive correlation of dose to epidemic severity and rates emerges in the GM model but not the AB model. Again, these differences can be attributed to the different methods of host contact (addressed above), which underscores the importance of host contact structure in quantitative predictions of caseloads. In contrast, epidemic rates are similarly negatively correlated to *P*(0) at 0 and moderate levels of *P*(0) (figure 4*c*,*d*), although the relative difference in rates is slightly more pronounced for the AB model relative to the GM model.

In both models, mean immunity in the population after the epidemic is lower when *P*(0) is higher, although again, the AB model predicts stronger effects. Interestingly, heterogeneity creates slightly higher population immunity relative to that resulting from a homogeneous population with *P*(0) equal to 1 (i.e. the mean of the heterogeneous case) (figure 5*a*,*b*). Further, note that a population with high partial immunity (*P*(0) = 2) has the same level of immunity as when the epidemic started, which is a result of the fact that in both models, the epidemic is very mild for all initial doses. The variance in immunity after an epidemic is generally lower in the GM model, but we do see agreement in the result that a population with heterogeneous *P*(0) will emerge from the epidemic with the largest amount of variation, a result that seems intuitive (figure 5*c*,*d*). However, less intuitive is the discrepancy between models in the levels of heterogeneity produced by different *P*(0) values: the AB model predicts that moderate *P*(0) will result in much higher levels of population variability in immunity relative to the levels emerging from a completely naive population, whereas the GM model predicts similarly low levels of variability in both cases. We attribute this effect to the contact structure of the AB model, which causes a significant proportion of a *P*(0) = 1 population never to become infected. These individuals who do not become infected increase the range of maximum *P*-values (figure 1*b*), resulting in larger variability in mean *P* than for the GM model (which sees no proportion of the *P*(0) = 1 population left uninfected). This increase in variability caused by a fraction of the population not being infected in the AB model is also evident in the heterogeneous *P*(0) case.

### (e) Control models

To compare the effect of dose-dependent replication and transmission against a system without such an effect, we created two models in which transmission was uncoupled from dose-dependent replication. In the constant-dose (CD) control model, the same dose is used to initiate all infections throughout the epidemic, whereas in the random-dose (RD) control model, the viral dose is random at every transmission, regardless of the viral load in the infected host. The control models give us a picture of the epidemic for two cases: one in which the within-host dynamics are constant and transmission is dose-independent, similar to a standard three-compartment SIR model, where dynamics of infected individuals depend on contact rates and infectious period, and the other in which the within-host dynamics are variable, as in our model, and transmission is dose-independent, which is similar to epidemiological models that include stochastic transmission (i.e. from random noise, not stochasticity from viral dose). Unlike the dose-dependent models, in the control models, the epidemic always progresses regardless of the value of *P*(0), since transmission does not involve a bottleneck (i.e. infected hosts always transmit to susceptible hosts) (electronic supplementary material, figures S1, and S3).

In the CD model, there is a slight decrease in epidemic severity and increase in epidemic rate with increasing dose (electronic supplementary material, figure S1). Note that an effect of dose emerges in this model because every infection is initiated with the same dose as the initial host received. Since a smaller dose results in a longer infection period, the maximum number of cases is larger and the time to reach this maximum is shorter. This effect is similar to that found in more standard epidemiological models that base transitions between compartments on infectious period. However, in the RD model, epidemic severity and rate are relatively insensitive to initial dose since the mean dose is constant across all initial dose conditions (electronic supplementary material, figure S3).

Across the control and dose-dependent models, we find that variation in immunity is highest for a heterogeneous *P*(0) population (electronic supplementary material, figures S2 and S4). Similar to the dose-dependent models, variation in immunity is increased for the AB model when not all susceptibles are infected in a *P*(0) = 1 or 2 population. The main difference between the control models and dose-dependent models is that variation in immunity changes with *P*(0) more in the dose-dependent model, emphasizing that dose-dependent replication and transmission enhances the effect of initial immunity in a population on immunity variation in the population after an epidemic.

There is also a significant qualitative difference between the control models and the dose-dependent model predictions of population immunity following the epidemic. In the dose-dependent model, the magnitude of effects of *P*(0) on epidemic behaviour is much higher than in the two control models (electronic supplementary material, figures S2 and S4). There is also a significant qualitative difference between the control models and the dose-dependent model predictions of mean population immunity following the epidemic. In the control models, mean immunity in the population is highest for moderate *P*(0) = 1, and then decreases for *P*(0) = 0 and *P*(0) = 2. On the other hand, in the dose-dependent model, we see that the order is *P*(0) = 0, *P*(0) = 1 and then *P*(0) = 2. Thus, when dose-dependent replication determines transmission, mean immunity in a population is highest in a completely naive population, whereas when transmission is uncoupled from dose-dependent within-host dynamics, mean immunity is highest for a population with a moderate initial immunity. Together, these differences indicate that dose-dependent within-host dynamics can affect host population immunity both quantitatively and qualitatively.

## 4. Discussion

We have compared two methods of nesting a dose-dependent within-host model of viral dynamics into an epidemiological model of a host population. The AB model provides a more complete integration of dynamics across these scales at the expense of computational intensity, whereas the GM model sacrifices some of the time-dependent quantities for computational ease. In analysing the effect of dose-dependent replication and transmission on population-level viral dynamics, we found that the GM model provides a reliable approximation of mean peak immunity in recovered hosts and viral doses throughout the epidemic time course. We also found that effects of initial immunity on the profile of population immunity after the epidemic are qualitatively different when dose-dependent replication and transmission dynamics are considered.

The GM model served as a mean-field approximation to the stochastic transmission in the AB model, leaving all common parameters equal. For this reason, the GM model is able to provide us similar insight to the behaviour of an epidemic in a population with varying levels of initial uniform immunity. Previous comparisons between deterministic compartmental models and AB models have shown that while AB models can provide a more realistic representation of epidemic dynamics because of their ability to capture more biological complexity and stochasticity, deterministic compartmental models can provide a reliable estimation of the mean behaviour (Brailsford & Hilton 2000; Schieritz & Milling 2003; Borschev & Filippov 2004; Morecroft & Robinson 2005; Demirel 2006; Rahmandad & Sterman 2008). This initial study revealed some method-based inconsistencies in predicted epidemic behaviour, e.g. epidemic rate and severity. Nevertheless, there are some strikingly congruent predictions across the models, and the GM model provides a promising foundation for the development of more complex analytical models that can examine long-term effects of specific immunity through reinfection of recovered individuals, as well as the evolution of immune escape.

Through comparison of the effect of initial population immunity on post-epidemic population immunity in the dose-dependent model and the two control models, we uncovered a qualitative disagreement between a system with dose-dependent transmission and systems with dose-independent transmission. The idea that levels of population immunity can be hindered by pre-existing immunity has been proposed for influenza epidemiology (Smith *et al*. 1999). This emphasizes the need for continued development of methods for understanding how within-host dynamics modify immunity in heterogeneous host populations, in order to inform disease control policies, such as identifying which groups should be given priority for vaccination. While several theoretical studies have addressed these issues using compartmental epidemiological models (Hethcote & VanArk 1987; Handeler & Muller 1996; Castillo-Chavez & Feng 1998; Hill & Longini 2003; Tennenbaum 2008), the impact of within-host dynamics on population immunity is yet to be considered in this context. Future work aimed at understanding how infection-induced variability in immunity levels impacts viral immune escape will help to understand how disease control methods can be implemented with the lowest risk of viral emergence.

In addition to within-host dynamics, host contact-network structure is a major determinant of between-host dynamics. In scale-free networks (i.e. the degree of connectivity between individuals is characterized by a power law distribution), epidemic spread is hierarchial, erupting first in the most connected individuals and last in the least connected ones (Barthélemy *et al*. 2004, 2005). Secondly, models that include this type of finite-network structure fail to exhibit epidemic threshold conditions that are typically observed in epidemiological models of finite populations (May & Lloyd 2001; Pastor-Satorras & Vespignani 2001). Future research should aim to integrate the effects of time-varying viral dynamics on transmission with those from network-constrained epidemic spread in order to predict how viral transmission is concurrently guided by within- and between-host factors. This type of integrative framework will reveal the constraints on viral emergence that can be used as a premise for developing effective methods of prevention.

Our calculation of *R*_{0} differed across the models. *R*_{0} represents the total number of new cases caused by a single infectious individual during an epidemic, and for this reason, the AB calculation of *R*_{0} is closer to that which could be calculated from real-world case data. However, the GM model does not allow for a calculation of this sort. Thus, as it stands, estimates of *R*_{0} using the GM model may not compare well with estimations from real-world data.

Most of the parameter values in our models were chosen based on observed data of acute viral infections, but some of the quantities are either not known or have been described only to a limited degree (e.g. quantitative details of the functional relationship of viral dynamics and specific immunity in primary versus secondary infections). Thus, we do not have an experimentally validated system, and this is an aspect of the research we hope can be improved in future work through collaboration with virologists and immunologists. Since we are mainly interested in model comparison in this work, however, we do not believe that this weakness affects our main finding that dose-dependent replication and transmission is an important driver of epidemiological dynamics, particularly the changing profile of host population immunity.

In future studies we plan to build on this framework by including mutation (i.e. multiple strains) to examine how these effects from dose-dependent within-host dynamics determine the emergence of immune-escape strains in a host population with changing immunity. Understanding the effects of dose-dependent within-host kinetics on mutation frequency and immune response is also important for developing viral disease therapies (e.g. drugs that impact viral kinetics) that optimize therapeutic gain with minimal risk of immune-escape evolution. For example, large doses would contain more viral mutants than smaller doses at each transmission event, which could result in more heterogeneity in specific immune responses within hosts in the population, but lower variation between hosts.

In conclusion, we have found that a mean-field approximation to discrete transmission events, the GM model, provides a reliable approximation for epidemic behaviour in terms of mean immunity. We have also found that dose-dependent replication and transmission have an effect on population immunity after an epidemic, in contrast to a dose-independent transmission model.

## Acknowledgements

This work was supported by NSF grant EF-0742373. K.M.P. was also supported by NIH grant R01 GM083983-01.

## Footnotes

One contribution of 14 to a Theme Issue ‘New experimental and theoretical approaches towards the understanding of the emergence of viral infections’.

- © 2010 The Royal Society