key: cord-0946844-rqqke40c authors: Heffernan, J.M.; Wahl, L.M. title: Improving estimates of the basic reproductive ratio: Using both the mean and the dispersal of transition times date: 2006-04-06 journal: Theor Popul Biol DOI: 10.1016/j.tpb.2006.03.003 sha: 78e6d72a5daf305f583b8fb8927eb222d9562e4d doc_id: 946844 cord_uid: rqqke40c In both within-host and epidemiological models of pathogen dynamics, the basic reproductive ratio, [Formula: see text] , is a powerful tool for gauging the risk associated with an emerging pathogen, or for estimating the magnitude of required control measures. Techniques for estimating [Formula: see text] , either from incidence data or in-host clinical measures, often rely on estimates of mean transition times, that is, the mean time before recovery, death or quarantine occurs. In many cases, however, either data or intuition may provide additional information about the dispersal of these transition times about the mean, even if the precise form of the underlying probability distribution remains unknown. For example, we may know that recovery typically occurs within a few days of the mean recovery time. In this paper we elucidate common situations in which [Formula: see text] is sensitive to the dispersal of transition times about their respective means. We then provide simple correction factors that may be applied to improve estimates of [Formula: see text] when not only the mean but also the standard deviation of transition times out of the infectious state can be estimated. The basic reproductive ratio, R 0 , is a fundamental concept in both within-host and epidemiological models of pathogen dynamics. Defined as the number of secondary infections caused by a single infectious individual in a susceptible population, R 0 is used to gauge the magnitude of the risk associated with a novel pathogen, and to estimate the degree of control that will be necessary to contain an outbreak. For example, in the last several years R 0 has been a key factor in our understanding of emerging diseases such as Severe Acute Respiratory Syndrome (SARS) (Choi and Pak, 2003; Lipsitch et al., 2003; Lloyd-Smith et al., 2003; Riley et al., 2003) and avian influenza (Stegeman et al., 2004; Ferguson et al., 2005; Longini et al., 2005) , livestock diseases such as Bovine Spongiform Encephalitus (Woolhouse and Anderson, 1997; Ferguson et al., 1999; de Koeijer et al., 2004) and Foot and Mouth Disease (Ferguson et al., 2001; Matthews et al., 2003) , as well as vector-borne disease such as Dengue (Luz et al., 2003) , Malaria (Hagmann et al., 2003; Smith and Ellis McKenzie, 2004) and West Nile Virus (Wonham et al., 2004; Cruz-Pacheco et al., 2005) . Recent analyses of HIV therapy (Smith, 2006) , river blindness (Filipe et al., 2005) and Escherichia coli in cattle (Matthews et al., 2006) also rely on this important concept (see Heffernan et al., 2005 for review) . Given precise knowledge of the probability distributions for times in the exposed or infectious states, the mathematical derivation of R 0 is straight-forward (Diekmann and Heesterbeek, 2000) . In practical terms, however, the estimation of R 0 from incidence or clinical data is far more difficult. Systems of ordinary differential equations (ODEs) or other simplifying arguments are used to derive estimates of R 0 in terms of measurable parameters (Dietz, 1975; Anderson and May, 1991; Dietz, 1993; Mollison, 1995; Diekmann and Heesterbeek, 2000; Hethcote, 2000; Brauer, 2002) . These estimates typically assume, however, that the lifetimes of individuals in each state are exponentially distributed. In an ODE formulation, this is a direct result of the Markov property that each state is memoryless, that is, the probability of making a transition to a new state does not depend on the length of time an individual has spent in the current state. For example, in all ODE models the probability of death is independent of the age of the individual or the time since infection; likewise the probability of recovery is independent of the time since infection. (For simplicity we use the word lifetime to mean the time spent in a single state or class such as ''infected''; we do not restrict this word to mean the entire lifetime of the individual.) As noted by many others, this assumption of constant transition probabilities is mathematically convenient, but is often difficult to defend on biological grounds (Lloyd, 2001a-c; Heesterbeek and Dietz, 1996; Perelson and Nelson, 1999; Nelson et al., 2004) . Previous studies have considered the effects of non-exponential lifetime distributions on the dynamics and persistence of childhood diseases (Grossman, 1980; Grenfell, 1997, 1998; Lloyd, 2001b,c) . Non-exponential distributions have also been used to describe the mean incubation period (Anderson, 1988; Lui et al., 1988) and infectious period (Blythe and Anderson, 1988; Castillo-Chavez et al., 1989; Mittler et al., 1998; Nelson et al., 2004) of Human Immunodeficiency Virus (HIV). For example, Lloyd (2001a) elucidated the sensitivity of R 0 to the probability distributions describing the viral life cycle. Unfortunately, for any practical study, the underlying probability distributions for transition times are almost certainly unknown. Nonetheless, in many cases the data at hand (or practical intuition) may provide more information than simply the mean lifetime. For example, while death from natural causes may be assumed to be exponential during an epidemic, recovery from the disease is almost certainly less dispersed. If the standard deviation of the recovery times is also known, we may be able to use this information to derive a more accurate estimate of R 0 . We might also improve our estimate of R 0 if the standard deviation is not known precisely, but is known to lie within a certain range. The goal of this paper is to derive improved estimates of R 0 for situations in which not only the mean, but some understanding of the dispersal of transition times is available to the clinical or epidemiological practitioner. To focus this investigation, we will restrict our attention to transition times out of the infectious state. We thus consider a large number of processes such as natural death, death due to disease, recovery and removal to quarantine, or for in-host models, loss of infected cells due to natural death, lysis, immune-mediated cell death, etc. The probability distributions of other parameters, such as infectivity, are also critical to disease dynamics (Nelson et al., 2004; Lloyd, 2001c) and we hope to address such issues in future work. We use the concept of ''competing processes'' to identify situations in which R 0 is sensitive to the dispersal of transition times. We then use the two most common lifetime distributions-the Weibull and the Gamma (see the Appendix)-to derive correction factors for R 0 , based on both the mean and the standard deviation of the appropriate lifetimes. The correction factor is modest ðp30%Þ when one of the competing processes is exponentially distributed, but approaches 200% when both processes are narrowly dispersed and have similar means. We motivate this study by first discussing several recent examples from the literature. 2. Recent estimates of R 0 are sensitive to lifetime distributions SARS: The estimation of R 0 for SARS received considerable attention Lipsitch et al., 2003; Lloyd-Smith et al., 2003; Riley et al., 2003; Nishiura et al., 2004) . Lipsitch et al. (2003) used an SEIR approach to derive R 0 ¼ b=ððm þ v þ wÞð1 À qÞÞ, where w is the mean daily rate at which infectious cases who are not in quarantine are detected and isolated, q is the fraction of persons contacted but not infected by an infectious person, b ¼ kb where k is the baseline daily number of contacts per capita and b is the probability of transmission per contact between a susceptible and an infectious, m is the per capita death rate due to disease and v is the per capita recovery rate. Fig. 1 of b such that R 0 ¼ 3, the mean value estimated by Lipsitch et al. (2003) , when lifetimes are distributed exponentially (dot). These estimates of R 0 are consistent with the range 2:2oR 0 o3:6 reported in the literature. We find that R 0 is quite sensitive to the dispersal of the relevant distributions. HIV epidemiology: R 0 for the epidemiological spread of HIV has also been derived from ODE models. For example, Hyman and Li (2000) derive R 0 ¼ rb=ðm þ vÞ, where b is the transmission probability per contact and r is the contact rate. The parameter v gives the rate of progressing to AIDS, while m gives the rate at which individuals are removed from the infectious compartment by other causes, for example, ceasing to be sexually active. The contact rate for individuals who have progressed to AIDS is assumed to be negligible; HIV is no longer transmitted at this stage. As an ODE formulation, this model assumes that the rate of progressing to AIDS is independent of the time since infection. In Fig. 2 (a) we demonstrate the sensitivity of R 0 to this assumption. Again we plot R 0 versus the CoV of a Weibull distribution of transition times to AIDS; the dot gives the case when these times are exponentially distributed. HIV immunology: R 0 is also central to the quantitative study of within-host HIV dynamics. The basic model of inhost infection consists of a set of ODEs that describe the interactions between susceptible cells, x, infected cells, y, and infectious virions v (Ho et al., 1995; Wei et al., 1995; Perelson et al., 1996; Nowak and May, 2000) : Susceptible cells are produced at a constant rate with mean l from a pool of precursor cells and die with rate d. Virions are budded with rate k and have probability Q of being infectious. Infectious virions can infect healthy cells via mass action dynamics where b describes the efficacy of this process. Infected cells die with rate a and virions are cleared at rate u. The basic reproductive ratio for this system depends on f ðtÞ, the probability that an infectious virion budded at time zero produces an infected cell that is still alive at time t. In terms of the underlying distributions we see that: whereûðrÞ denotes the probability that the virion is still in circulation at time r,ĝðrÞ denotes the rate at which a virion, which has existed for time r, infects andâðt À rÞ is the probability that an infected cell lives to be age t À r. Convolution is denoted by Ã. Substituting into Eq. (22) (see Appendix A.2) with bðtÞ ¼ kQ and assuming thatû,ĝ andâ are distributed exponentially, we find that R 0 ¼ kQbx 0 = ðaðbx 0 þ uÞÞ, where x 0 is the mean number of uninfected cells at the uninfected equilibrium. Eqs. (2) and (22) may also be used, however, when any of the underlying distributions are non-exponential. As an example, we consider the case when cell death time distributions are less dispersed than exponential, while all other processes are exponentially distributed. It is The dot denotes the case when these transition times are exponentially distributed; the distribution is described by a Weibull distribution with shape parameter r and mean 1=v (see text for details). All other lifetime distributions are assumed to be exponential. v ¼ 0:083 year À1 (Blower et al., 2001) , m ¼ 0:032 year À1 (Blower et al., 2001) and rb ¼ 0:23 year À1 individual À1 such that R 0 for the ODE model is 2 (Blower and McLean, 1994; Valesco-Hernandez et al., 2002) . (b) The initial growth rate r 0 for the in-host HIV model is plotted as a function of the CoV of the Weibull distribution of cell death times (dotted line); both uninfected and infected cells were assumed to have the same shape parameter although different means. Here r 0 is the solution to Eq. (23) when bðtÞ ¼ kQ and f ðtÞ ¼ R t r¼0 e Àðuþbx 0 Þr be ÀðsðtÀrÞÞ r dr, when s ¼ aGð1 þ 1=rÞ. The dots and error bars give r 0 AE one standard error of the mean as computed by Monte Carlo simulation (see . Parameters values: l ¼ 100 cells m l À1 day À1 , d ¼ 0:1 day À1 , b ¼ 0:0032 day À1 , a ¼ 0:5 day À1 , Q ¼ 0:05, k ¼ 62:5 day À1 and u ¼ 3 day À1 (Perelson and Nelson (1999) ; Wahl and Nowak (2000) ; Smith and Wahl (2004) ). straight-forward to demonstrate through Eq. (22) that the basic reproductive ratio is insensitive to this change. However, the initial growth rate of the virus strongly depends on these cell lifetime distributions, as demonstrated in Fig. 2(b) . This point will be taken up again in Section 4. Note that we have not considered any change in the infection time of an infectious virion since this case is not included in our restricted focus on distributions of transition times out of the infected state only. We likewise do not explore the delay between the infection of a cell and virion production (Herz et al., 1996; Culshaw and Ruan, 2000; Nelson and Perelson, 2002; Dixit and Perelson, 2004) . It is clear that changes to the mean times associated with transitions out of the infectious state(s) directly affect R 0 . The fact that R 0 is also sensitive to the shapes of these distributions, even if the mean is unchanged, is also relatively well-known (see Diekmann and Heesterbeek, 2000) . In the previous examples we found that R 0 sensitively depends on the shape of the distribution for the progression to AIDS. Likewise, for the SARS example, changes in the death time distribution affected R 0 . However, the shape of the death time distribution for infected cells had no effect on R 0 for in-host HIV infection. This can be explained using the concept of competing processes. For SARS, individuals leave the infected state either through recovery or death; recovery and death are competing, in this sense, to end the infectious period. In the epidemiological model of HIV, the infectious period ends if individuals progress to AIDS, or if they cease to be sexually active. In the in-host HIV model, however, the death of an infected cell is not in competition with another process. We define the ''effective infectious period'', L, as the time for which an individual remains infectious, irrespective of the process by which the individual leaves this compartment. When two (or more) processes are in competition to end the effective infectious period, not only the mean but the entire shape of the probability distribution for each process will affect the mean of L, and thus affect R 0 . Mathematically, consider two renewal processes which have transition times given by random variables A and B with means 1=a and 1=b, respectively. The probability distribution of A gives the distribution of times at which an infectious individual would leave the infectious compartment through process A; likewise for B. We assume that the two processes are independent, which is already implied in the system of differential equations. We also assume that the two processes are ''competing'', by which we mean that an infectious individual will leave the infectious compartment through one and only one of the two processes. In this case the survival probability, the probability that an individual is still in the infectious compartment at time t, can be written as since A and B are independent. Thus, f ðtÞ is simply the product of the tail distributions for A and B. The effective infectious period is then found by simply integrating f ðtÞ from zero to infinity. (far left), a ¼ b (middle) and bba (far right). When a ¼ b (closed circle), L is given by 1=ð2aÞ. When the means of the two distributions are very different (at the extremes to the left and right), L is dominated by one process and approaches either 1=a or 1=b. When a and b are of the same magnitude, however, the effective infectious period will be sensitive to both distributions. In the bottom panel we therefore examine the case when a ¼ b more closely. We now let the CoV of B vary, while A is exponentially distributed. B is described either by a Weibull distribution (solid line) or a Gamma distribution (dash dot line). When the CoV of B is 1, which corresponds to an exponential distribution, L ¼ 1=ð2aÞ (dot, middle inset). When the CoV is greater than 1, L decreases. When the CoV is less than 1 such that B is less dispersed than exponential (see left inset), L increases toward a plateau. For the Gamma distribution this asymptotic value of L is given by where k is the shape parameter, denoting the number of stages in the Gamma distribution. We also note that if aab (not shown), a similar limit is obtained, with L ! ð1=aÞð1 À e Àa=b Þ as the CoV of B narrows. For the Weibull distribution with shape parameter r, it is straight-forward to demonstrate that the same limits are approached as r ! 1. In fact as demonstrated in the figure, the Weibull and Gamma distributions give only slightly different effective infectious periods for the same CoV. We have found that the effective infectious period, and thus R 0 , is sensitive not only to distribution means but also to distribution shapes when two or more processes ''compete'' to end the infectious period, and when the mean times for the competing processes are of the same order of magnitude. As discussed previously, methods for estimating R 0 from incidence or clinical data typically assume that all the underlying processes are exponentially distributed. In the following section, we derive more accurate estimates for R 0 , assuming that one or both of the competing processes are non-exponentially distributed. Rather than re-deriving R 0 for a number of models (SIR, SEIR, vector etc.), we report a ''correction factor'', f: the ratio of R 0 when the lifetimes are non-exponentially distributed to the value R 0 that would be calculated assuming exponential lifetimes. Where possible we report f as a function of the CoV of the relevant distribution; in many cases we are able to derive limiting values of f as this distribution narrows. We also find that f can approach this asymptotic limit very quickly, and simple formulae can be used to correct estimated values of R 0 , or to gauge the sensitivity of R 0 to dispersion in the underlying distributions. Throughout the following section we again assume that two renewal processes compete to end the infectious period. The transition times for these processes are given by random variables A and B, with means 1=a and 1=b, respectively. As stated previously, we assume throughout this derivation that other processes involved in the infection time course (for example, entering the infectious state) occur at constant rates. This implies, in particular, that infectiousness is constant throughout the infectious period, such that R 0 may be written as R 0 ¼ bL, where L is the effective infectious period and the transmission rate b is constant. We are then able to estimate f as the ratio of the true value of effective infectious period, L, to the value of L predicted by an ODE model, 1=ða þ bÞ. In this section we assume that only random variable A is exponentially distributed. For example, it may be natural to assume that the recovery distribution is less dispersed than exponential, while the natural death rate is roughly constant over the time course of the infection. (In the latter case the correction factor would be modest, however, because the means of the two distributions are usually quite different.) Removing infectious individuals to quarantine is another process which competes with recovery and is more likely to have a mean of the same order of magnitude. Assume that A is exponentially distributed and B is either described by a Gamma distribution or a Weibull distribution. If B is described by a Gamma distribution with k stages, the effective infectious period can be found by integrating Eq. (3): and thus the correction factor is If a ¼ b, which can be used to approximate the case when a and b are similar, this factor reduces to Of course in practice we are unlikely to have in hand an estimate of k. Recall, however, that the CoV of the Gamma distribution is given by . Thus if C can be estimated, we can write the correction factor as It is also quite useful to note that, in analogy with Eq. (4), which reduces to f ¼ 2ð1 À e À1 Þ % 1:26 if a ¼ b. In practice this limit is approached very quickly as the distribution of B becomes less dispersed, as shown in Fig. 4 . As a practical example, assume as in an ODE model that the rate at which individuals are moved to quarantine is independent of the time since infection, with a mean time of about 1 week. If we know in addition that the recovery time is 1 week AE2 days, the CoV is about 0.3, and f is within about 2% of the limit given in Eq. (9). The approximation is even better if the CoV is smaller than 30%. Fig. 4 plots the correction factor for different ratios of the distribution means, assuming B is described by a Gamma distribution. We see that f decreases with the CoV of B for all values of b=a. Note that f is maximized at % 1:3 when b=a ¼ 1:79. This can be found numerically by solving e a=b ¼ 1 þ ða=bÞ þ ða=bÞ 2 , which is obtained by differentiating Eq. (9). If B is assumed to be a Weibull distribution with shape parameter r we obtain: where GðxÞ is the Gamma function. This can only be solved numerically, however results (not shown) are extremely close to those illustrated in Fig. 4 for the Gamma distribution. In many cases, it may be more realistic to assume that neither of the processes in competition are distributed exponentially. For example, both recovery and death due to disease typically depend on the time that a host has been infected. The correction factor for such cases is more complex. If A and B are described by Weibull distributions with shape parameters r a and r b , respectively, we find Similarly, when A and B are Gamma distributions with shape parameters k a and k b we find Both of these expressions must be determined numerically. These expressions can be simplified considerably, however, if the shapes of the two distributions are known to be similar. Recall, for example, that a shape parameter r ¼ 3 in the Weibull distribution corresponds to a distribution that is roughly Gaussian, with a CoV of around 30-40%. If processes A and B can be approximated by Weibull distributions with the same shape parameter r, we find This allows for three further simplifications: when, in addition to the shapes, the means of the two processes are roughly equal ða % bÞ we find f % 2 ðrÀ1Þ=r . If the means differ but the CoV is small (say less than about 30%), we can take the limit as r ! 1 to find f % ða þ bÞ=maxða; bÞ. Finally, if the means are similar and the CoV is small, we have f % 2. These results are summarized in Table 1 . Although it may seem unlikely that the two relevant distributions would have exactly the same shape, the assumption that both shapes have a CoV of less than 30% may be met in many cases. Similarly, if A and B can be approximated by Gamma distributions with different means but the same shape parameter k, the correction factor is where Hða; b; c; xÞ is the generalized hypergeometric function. Although this expression is still unwieldy, when a % b this reduces nicely to since Hð½Àk; 1; k; À1Þ ¼ 0. Once again although k may not be known, we can substitute k ¼ 1=C 2 , where C is the CoV of the distribution, to find f as a function of the CoV. It is assumed that A is exponentially distributed with mean 1=a. Fig. 5 plots the correction factor, determined numerically, versus the CoV of A and B when A and B are both described by Weibull distributions but do not have the same shape. In the top panel the distribution means are the same, b ¼ a; in the middle panel b4a. In each case the correction factor decreases as the CoVs of A and B increase. Note once again that the limit for small CoV gives a good approximation for f whenever the CoV is less than about 30%. In the bottom panel we plot the correction factor for the analytically simpler case when the distributions have similar shapes. We take the ratio of the means, b=a to be between 1 and 5. Again we see that f is very close to its limiting value ða þ bÞ=maxða; bÞ whenever the CoV is small. The larger b=a is, the more quickly this very simple limit is approached. Results for the Gamma distribution (not shown) were once again almost indistinguishable. Our main results are summarized in Table 1 . We find that simple closed-form expressions for the correction factor are available when one of two competing processes is described by an exponential distribution, while the other has a measurable CoV or can be assumed to have a CoV less than about 30% (first row of Table 1 ). If neither competing process is distributed exponentially, simple approximations are only possible when the two processes have roughly the same shape, although their means may of course differ (second row). The bold entries in the table give even simpler results for the case when the two competing processes have similar means. In the limiting case when the processes have similar means and similar shapes, we see that the correction factor approaches 2 as the shapes of the distributions narrow. These derivations assume that any non-exponential distributions can be approximated by either a Gamma or a Weibull distribution, and that some knowledge of the dispersal of the distribution about its mean is available. We found very little difference in the effective infectious period when either a Gamma or Weibull was assumed (see Fig. 3 ). Interestingly, each of these distributions had analytical advantages, and reduced to simple closed-form expressions in different cases. We thus include results from both distributions in Table 1 where appropriate. Table 1 Approximations to the correction factor, f, for competing processes A and B with means 1=a and 1=b, respectively A, B non-exponential ða þ bÞða r þ b r Þ À1=r a þ b maxða; bÞ same shape 2 ðqÀ1Þ =q, or 2 1 À 2 1À2k Gð2kÞ kGðkÞ 2 2 C gives the CoV of the non-exponential distribution; for the Gamma distribution k ¼ 1=C 2 . Entries in bold are simplifications for the case a ¼ b. The formulae we derive are applicable to situations in which two processes compete to end the infectious period. Examples in epidemiology include death due to disease, recovery, natural death and quarantine; for in-host models, examples include the loss of infected cells due to natural death, lysis, or immune-mediated cell death. These formulae apply both to models with direct transmission from one host to another e.g. SARS, or to vector models where the pathogen is carried between hosts e.g. malaria. It is also straightforward to apply the correction factors provided here to situations in which more than two processes are in competition (as in the SARS example, Fig. 1) . The dispersal of lifetimes involved in entering the infectious state may also affect R 0 . For an SEIR model, a slightly different correction factor can be easily derived for transitions into the infected state from the exposed state. Likewise for a virus or other types of in-host pathogens, the transition from one infected state to another may involve transitions both out of and into an infected state. For example, in the in-host model of HIV, the infection of a cell is a transition out of the infectious state for a virion, but it also represents a transition into the infected state for an infected cell. These situations suggest clear avenues for future work. In future we also hope to relax the assumption that the transmissibility of the disease is constant throughout the infectious period. We would also like to discuss, briefly, the effects of lifetime distributions for ''non-competing'' processes. For example, in Section 2 we found that lifetime distribution for the infected cells did not affect R 0 . This is because the death of infected cells is not in competition with any other process. Nonetheless, solving system (1) for the infected equilibrium we find where x 0 is the initial number of uninfected T-cells (Nowak and May, 2000; Lloyd, 2001a) . If, instead, we assume that the death time distributions for the uninfected and infected cells are distributed as arbitrary distributionsdðtÞ andâðtÞ with means 1=d and 1=a, respectively, the infected equilibrium is given bȳ Thus, we find that the infected equilibrium will change even though R 0 is not affected when the death of infected cells is non-exponentially distributed. We also find that ifx and x 0 are known, R 0 can be directly calculated from disease data irrespective of this distribution. If the distributions of competing processes are non-exponential, however, the relationship between R 0 and the infected equilibrium will change. For non-competing processes, the initial growth rate r 0 may be sensitive to lifetime distributions even though R 0 is unaffected, as shown in Fig. 2. Fig. 6 plots the time distribution of infectious virions budded from a single infected cell (bðtÞf ðtÞ from Eq. (22)) when the distribution of cell lifetimes is exponential (dashed line) or approximately Gaussian (solid line). In both cases the area under the curve is the same, and thus R 0 does not change. When lifetimes are Gaussian-distributed, however, the infected cells bud the same number of virions in a shorter period of time, thus increasing the initial growth rate (Fig. 2(b) ). In this case, the conventional set of differential equations (1) allows us to obtain the correct value for the basic reproductive ratio, but does not model the initial growth rate or subsequent dynamics accurately. The area under each curve is the same i.e. R 0 is the same, but r 0 changes. Parameters values are the same as in Fig. 2. Like the infected equilibrium, r 0 can be determined from disease data for many infectious pathogens, and this value can be used to determine R 0 (see Hethcote and Tudor, 1980; Nowak and May, 2000; Lipsitch et al., 2003; Lloyd, 2001a-c for examples) . We note that this method for estimating R 0 is not robust in the presence of nonexponential lifetime distributions for either competing or non-competing processes. This work was supported by the Natural Sciences and Engineering Research Council of Canada. We are indebted to Baolai Ge and Gary Molenkamp for their invaluable technical expertise. The Weibull distribution, first introduced by Wallodi Weibull to model failure times in reliability applications (Weibull, 1951) , is one of the most widely used distributions in survival analysis. The Weibull distribution is also a natural model for lifetime (survival) probabilities in many instances, as it models waiting times until a critical event occurs. We employ the two parameter Weibull distribution: where r is the shape parameter of the distribution and s ¼ ð1=mÞGð1 þ 1=rÞ is the scale parameter. Here GðxÞ is the gamma function and m is the mean of the distribution. Intuitively, this gives the distribution of times until a critical event occurs, such as cell death, the clearance of an infection or the transition from exposed to infectious. If r41 the failure rate is increasing, that is, the probability that an individual will leave the state increases over time; if ro1 the failure rate is decreasing and if r ¼ 1 the failure rate is constant. A shape parameter of one corresponds to an exponential distribution and a shape parameter of three corresponds to a distribution that is approximately Gaussian. The CoV for the Weibull distribution with shape parameter r is ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Gð1 þ 2=rÞ=Gð1 þ 1=rÞ 2 À 1 q . The Gamma distribution is also used to model failure times in reliability applications. We use the two parameter Gamma distribution given by the following: where l À1 ¼ m=k is a scale parameter, 1=m is the mean and k is the shape parameter. This distribution, like the Weibull, includes the exponential distribution when k ¼ 1. The Gamma distribution is not as widely used as the Weibull distribution in reliability analysis, but it is commonly used to describe stages of infection in epidemic models since it is simply a sum of exponential distributions; the sums of independent and identically distributed exponential random variables have a Gamma distribution. The Gamma distribution, unlike the Weibull, cannot be used to approximate a normal distribution. The CoV for the Gamma distribution with shape parameter k is equal to 1= ffiffiffi k p . A.2. The basic reproduction ratio, R 0 Consider a large population and let f ðaÞ be the probability that a newly infected individual remains infectious for at least time a. This is called the survival probability. Also, let bðaÞ denote the average number of newly infected individuals that an infectious individual will produce per unit time when infected for total time a. Then, R 0 is given by See Heesterbeek and Dietz (1996) for a detailed description of this equation and a historical overview. The definition of f ðaÞ must be extended for cases in which a series of states are involved in the ''reproduction'' of an infected individual. For example, in within-host viral dynamics, an infectious virion produces an infected cell which produces more infectious virions; this complete cycle must be taken into account in our derivation of R 0 . One approach is to define f ðaÞ in this case as f ðaÞ ¼ Z a 0 probðinfectious virion produced at time 0 exists at time tÞ Á probðinfectious virion of age t infects a cellÞ Á probðinfected cell lives to be age a À tÞ dt, and bðaÞ is simply the average number of infectious virions produced by a cell which has been infected for time a. Note that this extension to the R 0 concept is typically assumed in epidemiology and HIV dynamics (see , for review), but differs from the standard Next Generation Approach (Diekmann et al., 1990; Diekmann and Heesterbeek, 2000; van den Driessche and Watmough, 2002) . A.3. The initial growth rate, r 0 While the basic reproductive ratio gives the overall number of new infections per infected individual, it gives no measure of when these infections occur. In contrast r 0 gives the number of new infections per infected individual per unit time, in a susceptible population (Heesterbeek and Dietz, 1996) . Using the survival probability f ðaÞ, as derived above, we propose that r 0 can be determined by solving the following renewal equation: Here yðtÞ denotes the incoming rate of new infected individuals, or the incoming rate of individuals to state 1 of an infection cycle, at time t. Eq. (23) thus implies that new infections at time t are produced by new infections at time t À s which are still infectious s time units later. The dðtÞ term is a Dirac delta function, and is necessary because we wish to determine r 0 after the instantaneous input of one infectious individual at time 0. Conveniently, Eq. (23) is in the form of a convolution integral, which lends itself to solution by a Laplace transform. Taking the Laplace transform of both sides we obtain Y ðoÞ ¼ Y ðoÞF ðoÞ þ 1, where F ðoÞ and Y ðoÞ denote the Laplace transforms of bðsÞf ðsÞ and yðtÞ, respectively. Isolating Y ðoÞ and taking the inverse Laplace transform we find that Thus, the initial growth rate r 0 can be found by determining the largest exponent in the inverse Laplace transform for yðtÞ. The epidemiology of HIV infection: variable incubation plus infectious periods and heterogeneity in sexual activity Infectious Diseases of Humans: Dynamics and Control Prophylactic vaccines, risk behaviour change and the probability of eradicating HIV in San Francisco Live attenuated HIV vaccines: predicting the tradeoff between efficacy and safety Distributed incubation and infectious periods in a model of the transmission dynamics of the human immunodeficiency virus (HIV) Basic ideas of mathematical epidemiology On the role of long incubation periods in the dynamics of acquired immunodeficiency syndrome (AIDS) A simple approximate mathematical model to predict the number of severe acute respiratory syndrome cases and deaths Modelling the dynamics of West Nile Virus A delay-differential equation model of HIV infection of CD4 T-cells Quantifying BSE control by calculating the basic reproduction ratio R 0 for the infection among cattle Mathematical Epidemiology of Infectious Diseases: model Building, Analysis and Interpretation On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases Transmission and control of arbovirus diseases The estimation of the basic reproduction number for infectious diseases Complex patterns of viral load decay under antiretroviral therapy: influence of pharmacokinetics and intracellular delay Epidemiological determinants of spread of causal agent of severe acute respiratory syndrome in Hong Kong Estimation of the basic reproduction number of BSE: the intensity of transmission in British cattle The foot and mouth epidemic in Great Britain: pattern of spread and impact of interventions Strategies for containing an emerging influenza pandemic in Southeast Asia Human infection patterns and heterogeneous exposure in river blindness Oscillatory phenomena in a model of infectious diseases. Theor Malaria and its possible control on the Island of Principe The concept of R 0 in epidemic theory Monte Carlo estimates of natural variation in HIV infection Perspectives on the basic reproductive ratio The mathematics of infectious diseases Integral equation models for endemic infectious diseases Viral dynamics in vivo: limitations on estimates of intracellular delay and virus decay Rapid turnover of plasma virions and CD4 lymphocytes in HIV-1 infection An intuitive formulation for the reproductive number for the spread of diseases in heterogenous populations Disease extinction and community size: modeling the persistence of measles Effect of variability in infection period on the persistence and spatial spread of infectious diseases Transmission dynamics and control of sever acute respiratory syndrome The dependence of viral parameter estimates on the assumed viral life cycle: limitations of studies of viral load data Destabilization of epidemic models with the inclusion of realistic distributions of infectious periods Realistic distributions of infectious periods in epidemic models: changing patterns of persistence and dynamics. Theor Curtailing transmission of severe acute respiratory syndrome within a community and its hospital Containing pandemic influenza at the source A model-based estimate of the mean incubation period for AIDS in homosexual men Uncertainties regarding dengue modeling in Rio de Janeiro Neighbourhood control policies and the spread of infectious diseases Heterogeneous shedding of Escherichia coli O157 in cattle and its implications for control Influence of delayed viral production on viral dynamics in HIV-1 infected patients Epidemic Models: Their Structure and Relation to Data Mathematical analysis of delay differential equation models of HIV-1 infection An age-structured model of HIV infection that allows for variations in the production rate of viral particles and the death rate of productively infected cells Modelling potential responses to severe acute respiratory syndrome in Japan: the role of initial attack size, precaution, and quarantine Virus Dynamics Mathematical analysis of HIV-1 dynamics in vivo HIV-1 dynamics in vivo: virion clearance rate, infected cell lifespan, and viral generation time Transmission dynamics of the etiological agent SARS in Hong Kong: impact of public health interventions Statistics and dynamics of malaria infection in Anopheles mosquitoes Adherence to antiretroviral HIV drugs: how many doses can you miss before resistance emerges? Distinct effect of protease and reverse transcriptase inhibition in an immunological model of HIV-1 infections with impulsive drug effects Avian influenza A virus (H7N7) epidemic in The Netherlands in 2003: course of the epidemic and effectiveness of control measures Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission Adherence and drug resistance: predictions for therapy outcome Viral dynamics in HIV-1 infection A statistical distribution function of wide applicability An epidemiological model for West Nile virus: invasion analysis and control applications Understanding the epidemiology of BSE