Journal of Theoretical Biology 395 (2016) 97–102
Contents lists available at ScienceDirect
Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/yjtbi
A theoretical framework to identify invariant thresholds in infectious disease epidemiology M. Gabriela M. Gomes a,b,c,n, Erida Gjini d, Joao S. Lopes d, Caetano Souto-Maior d, Carlota Rebelo e a
CIBIO-InBIO, Centro de Investigação em Biodiversidade e Recursos Genéticos, Universidade do Porto, Portugal Instituto de Matemática e Estatística, Universidade de São Paulo, Brazil c Liverpool School of Tropical Medicine, United Kingdom d Instituto Gulbenkian de Ciência, Oeiras, Portugal e Departamento de Matemática e Centro de Matemática, Aplicações Fundamentais e Investigação Operacional, Universidade de Lisboa, Portugal b
H I G H L I G H T S
Host heterogeneity modifies the relationship between R0 and disease prevalence. Neglecting heterogeneity results in underestimated efforts to meet control targets. Invariant transmission thresholds are robust indicators for control planning.
art ic l e i nf o
a b s t r a c t
Article history: Received 14 April 2015 Received in revised form 17 January 2016 Accepted 20 January 2016 Available online 8 February 2016
Setting global strategies and targets for disease prevention and control often involves mathematical models. Model structure is typically subject to intense scrutiny, such as confrontation with empirical data and alternative formulations, while a less frequently challenged aspect is the widely adopted reduction of parameters to their average values. Focusing on endemic diseases, we use a general transmission model to explain how mean field approximations decrease the estimated R0 from prevalence data, while threshold phenomena – such as the epidemic and reinfection thresholds – remain invariant. This results in an underestimation of the effort required to control disease, which may be particularly severe when the approximation inappropriately places transmission estimates below important thresholds. These concepts are widely applicable across endemic pathogen systems. & 2016 Elsevier Ltd. All rights reserved.
Keywords: Endemic infection Epidemic threshold Reinfection threshold Heterogeneity Global health
1. Introduction Traditional epidemiologic research classifies individuals by similarity of stipulated characteristics and conceives models for disease distribution in a population that is compartmentalized (Rothman, 2012). In infectious diseases, hosts move between compartments according to processes, such as disease progression and transmission, describing dynamic patterns of disease in the population (Anderson and May, 1991). The logic of compartmental study design motivates compartmental models calibrated by average parameters: the “mean field” approximation. Of particular interest from such models are threshold parameters, such as the n Corresponding author at: CIBIO-InBIO, Centro de Investigação em Biodiversidade e Recursos Genéticos, Universidade do Porto, Portugal. E-mail address: [email protected] (M.G.M. Gomes).
http://dx.doi.org/10.1016/j.jtbi.2016.01.029 0022-5193/& 2016 Elsevier Ltd. All rights reserved.
basic reproduction number R0, and the effective reproduction number Re, which have become instrumental in the design and evaluation of control strategies (Diekmann et al., 1990; Becker et al., 1995; Anderson and May, 1991). Threshold parameters, generally dependent on model structure, can only be estimated indirectly and subject to assumptions whose impact is often unrecognized (Heffernan et al., 2005). Here we assess the behavior of transmission threshold parameters under heterogeneity in individual capacity to acquire immunity upon infection and in the rate at which contacts are made with other individuals. Populations are undoubtedly composed of individuals that differ in their propensity to acquire infections and in their potential to transmit to others. Overall, variation results from a mix of genetic and environmental factors, including social and physical aspects. Due to demographic processes, differential selection is likely to occur within groups, resulting in patterns of disease for the population as a whole that differ from the
98
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
expectation (Proscham and Sethuraman, 1976; Ball, 1985; Vaupel and Yashin, 1985). For this reason, variation in individual characteristics and group structure must be considered throughout the entire research process that establishes immunization practices, as recognized almost 50 years ago in a classical paper by Fox et al. (1971). Today there is an extensive literature dealing with mathematical analyses of how herd effects of vaccination programs help to reduce incidence and prevalence of infections among the human population. Until recently, most of this work has assumed that within each compartment, individuals are equally susceptible to infection, equally infectious if infected, and equally active in their social contacts. Although deviations of heterogeneous systems from mean field approximations have been well characterized mathematically for epidemic scenarios (Boylan, 1991; Andreasen, 2011; Clancy and Pearce, 2013; Katriel, 2012; Novozhilov, 2012), translations and applications in epidemiology research, especially with regards to endemic pathogens, remain limited. In this paper, we recall a few exceptions that are worth retaining before proceeding with general considerations on the relationship between threshold parameters and disease prevalence. Implication for disease control and elimination are discussed in Section 4.
studies implicate the droplet nucleus mechanism in the transmission of tuberculosis, measles, influenza and many other infections of the respiratory tract (Riley, 1974). Over the years several authors have independently challenged the mean field approximation in TB transmission systems. Murphy et al. (2002) evoked different genetic susceptibility distributions, demographic factors, and transmission intensity, to explain the wide variation in endemic TB levels between countries. Gomes et al. (2012) analyzed published data relating the reinfection proportion among recrudescences and TB incidence to estimate the decomposition into low and high risk groups that could best reproduce the worldwide dataset. Similarly, an independent study based in the city of Rio de Janeiro, Brazil (Dowdy et al., 2012), estimated relative contributions of high and low risk groups compatible with global estimates. Heterogeneity in propensity of infectious individuals to cause secondary cases has also been postulated to explain the observed skewness in genotypic cluster size distributions of M. tuberculosis (Vynnycky et al., 2001; Ypma et al., 2013) and outbreak dynamics in a number of directly transmitted diseases (Lloyd-Smith et al., 2005).
3. Endemic diseases and threshold parameters 2. Heterogeneity and transmission routes 2.1. Sexual transmission Since the beginning of AIDS epidemic and the identification of the aetiological agent in the 1980s, it has been evident that heterogeneity in individual sexual behaviors must be considered in the interpretation of epidemiological data (Anderson et al., 1986; Hyman and Stanley, 1988). The notion is so consensual that most textbooks on infectious disease epidemiology written since then have a major chapter or section dedicated to heterogeneity in sexually transmitted diseases. Much research has been devoted to estimating such contact networks in diverse settings and with different methods in order to accurately predict HIV transmission dynamics (Leigh Brown et al., 2011; Jones and Handcock, 2003). 2.2. Vector transmission Relatively recently, Smith et al. (2005) developed a mathematical framework that allows the estimation of heterogeneity in malaria infection rates from the relationship between parasite rates and entomological inoculation rates in multiple populations. The model was fitted to data from more than 90 localities, a distribution of individual susceptibilities was estimated, and the implications for disease control were discussed for this vectortransmitted disease. Notably, this heterogeneity was later found to greatly enlarge the range of R0 for malaria transmission (Smith et al., 2007), and entered the discussion of control strategies aiming at R0 o 1. 2.3. Airborne transmission The study of heterogeneity in airborne diseases has been less explored since parameter estimation becomes more difficult, although heterogeneity in social contacts, for example, has started to be recognized (Mossong et al., 2008). A prominent case in the topic is tuberculosis (TB), a disease transmitted by droplet nuclei carrying Mycobacterium tuberculosis secreted from patients with pulmonary TB when coughing, sneezing, spitting and other respiratory acts. These tiny droplets disperse throughout the air of enclosed spaces, such as rooms and buildings, and transmission occurs when they are inhaled by another person. Epidemiologic
The approaches used by Smith et al. (2005) in malaria, and Gomes et al. (2012) in tuberculosis, have a similar strategy of fitting a relatively simple model simultaneously to epidemiological data from multiple populations to infer common underlying heterogeneities. This enables the simultaneous estimation of key threshold parameters and the exploration of implications for disease control in specific settings (Gomes et al., 2014). We proceed by defining two threshold parameters of interest in the context of the mean field susceptible-infected-recovered (SIR) model, and then by analyzing the effects of adding heterogeneity incrementally. 3.1. Mean field SIR To define the threshold parameters and their behavior under individual heterogeneity, we rely on a minimal SIR model, which describes the rates of change in the proportions of the population that are susceptible (S), currently infected (I), and recovered (R) by the system of differential equations: dS ¼ μ λS μ S dt dI ¼ λðS þ σ RÞ γ þ μ I dt dR ¼ γ I λσ R μR; dt
ð1Þ
where λ ¼ βI is the per capita rate of infection (force of infection) given an effective contact rate, β, between susceptible and infected individuals, γ is the rate at which individuals recover from infection, μ is the birth rate, here assumed equal to the death rate, and σ is the relative susceptibility of recovered individuals with respect to those who have never encountered the pathogen. The system has two threshold parameters: 1) The epidemic threshold (Becker et al., 1995), above which epidemics occur and infection remains endemic, is given by R0 ¼ 1, where R0, defined as the average number of secondary infections generated by one infectious individual in a totally susceptible population, is given by R0 ¼ β = γ þ μ ; 2) The reinfection threshold (Gomes et al., 2004, 2005), above which epidemics occur and infection remains endemic in a population where all individuals have previously experienced
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
99
0
0.6
−2
10 0.4 0.2
−4
0
0
10
10
1
10
0
10
proportion infected
0.06
0.04 −2
10 0.02
0
proportion infected
proportion infected
0.8
proportion infected
10
1
−4
0
1
2
3
4
5
0
basic reproduction number, R
0
0.5
1
susceptibility
0
10
10
1
10 basic reproduction number, R
0
Fig. 1. Endemic equilibrium as a function of R0 and contact heterogeneity. Each colored curve assumes a different distribution of acquired immunity as represented on the right, with the same mean susceptibility hsi ¼ 0:5, but different variance. A, B Contact rates are distributed uniformly [model (2)]. C, D Contact rate distribution assumes some individuals have more contacts (higher risk) than others [model (7) with p1 ¼ 0:96, κ 1 ¼ 0:0875 and p2 ¼ 0:04, κ 2 ¼ 22:9, supporting a variance-to-mean ratio of 20]. Vertical lines mark the transmission thresholds: epidemic threshold R0 ¼ 1 (solid); reinfection threshold Re ¼ ⟨s⟩R0 ¼ 1 (dashed). Colored dashed lines represent the cases where all individuals have previously experienced the pathogen and become sensitized. Left and right panels represent the same curves in linear and logarithmic scales, respectively. The remaining parameters were set at: γ¼ 1 yr 1 and μ¼ 1/70 yr 1.
the pathogen and become partially immunized, is Re ¼ 1, where Re is the effective reproduction number given by Re ¼ σ R0 . The endemic infection prevalence at the equilibrium of the mean field SIR model is represented in terms of R0 by the solid red line in Fig. 1A,B (in linear and logarithmic scales, respectively), and the epidemic and reinfection thresholds are represented by the vertical solid and dashed lines, respectively. 3.2. SIR with distributed acquired immunity To assess the behavior of these thresholds under heterogeneity in the level of immunity acquired by exposure to the pathogen, we formulate an extended model in which the relative susceptibility of individuals in the R compartment, denoted by s and defined between 0 and 1, is described by a probability density function, R1 qðsÞ with mean ⟨s⟩ ¼ 0 qðsÞsds. In this sub-section, which refers to Fig. 1A,B and to the Supplementary material, Section S1, we analyse the model extension just described, formally written as: dS ¼ μ λS μS dt dIðsÞ ¼ λðqðsÞS þ sRðsÞÞ γ þ μ IðsÞ dt dRðsÞ ¼ γ IðsÞ λsRðsÞ μRðsÞ; ð2Þ dt R 1 where λ ¼ β 0 IðsÞds is the force of infection upon naive individuals. The basic reproduction number for this model has the same dependence on basic parameters as the mean field approximation, R0 ¼ β = γ þ μ , and so does the epidemic threshold, since the reinfection process, by definition, does not affect R0. As for the effective reproduction number, it is now written as Re ¼ hsiR0 , showing that the expression for the reinfection threshold is also unaffected by heterogeneity in reinfection since only the mean susceptibility of recovered individuals appears in the formula. We represent endemic equilibrium infection prevalence in terms of R0, and compare the described trends assuming that qðsÞ is a beta
distribution with mean hsi ¼ σ . We choose a beta distribution because it fits the criterion of being defined in the interval (0,1) and allows a spectrum of shapes ranging from homogeneous to all-or-nothing (central panels in Fig. 1). In addition, we consider two limiting discrete scenarios: a homogeneous implementation where all previously infected individuals recover with the same factor of partial protection, qðsÞ ¼ δðs σ Þ; for heterogeneous we adopt a polarized (all-or-nothing) implementation where some recovered individuals retain total protection after infection while the remaining retain none, qðsÞ ¼ ð1 σ ÞδðsÞ þ σδðs 1Þ. The results using beta distributions are obtained numerically, whereas the discrete limits are treated analytically as described below. The model with homogeneous acquired protection is written as (1) above. This previously analysed system (Gomes et al., 2004, for example) admits a disease-free equilibrium for any value of the parameters. This state is stable when R0 o 1, and undergoes a transcritical bifurcation at R0 ¼ 1, losing its stability and giving rise to a stable endemic equilibrium when R0 4 1. The endemic equilibrium for system (1) can be solved explicitly, and in particular the proportion infectious is given by: I¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 ffi hsi R0 γ þ μ μ γ þ μ þ hsi2 R0 γ þ μ μ þ γ þ μ 2hsi R0 γ μ þ μ γ þ μ : 2hsiR0 γ þ μ
ð3Þ The model with all-or-nothing acquired protection is written as: dS ¼ μ λS μ S dt dIð0Þ ¼ ð1 σ ÞλS γ þ μ Ið0Þ dt dIð1Þ ¼ λðσ S þ Rð1ÞÞ γ þ μ Ið1Þ dt dRð0Þ ¼ γ Ið0Þ μRð0Þ dt dRð1Þ ¼ γ Ið1Þ λRð1Þ μRð1Þ dt
ð4Þ
where λ ¼ βI is the force of infection for I ¼ Ið0Þ þ Ið1Þ, and σ is the
100
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
proportion of individuals who fail to retain protection to subsequent infections (susceptibility factor 1) while 1 σ is the proportion acquiring total protection (susceptibility factor 0). Solving system (4) for the endemic equilibrium we find that the proportion infectious is given by: I¼
hsiR0 γ þ ðR0 1Þμ
γ þμ þ
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 hsiR0 γ þ ðR0 1Þμ þ γ þ μ 4hsiR0 γ γ þ μ : 2R0 γ þ μ
ð5Þ The solutions given by expressions (3) and (5) are represented in Fig. 1A,B by the solid red and blue lines, respectively, as functions of R0 . Homogeneous protection maintains prevalence levels that are consistently higher than those maintained under the heterogeneous version with the same mean. Two intermediate scenarios, with continuously distributed acquired protection (system (2)), were treated numerically by discretizing the beta distributions and solving the resulting ordinary differential equations in MATLAB (yellow and green in Fig. 1A,B). Noticeably, partial protection upon recovery generates the same reinfection threshold at R0 ¼ 1=hsi in all cases. Following Gomes et al. (2005), further insight into the reinfection threshold is obtained by considering a sub-model of system (2), in which the S compartment is bypassed and all individuals enter the population into the respective R classes as if they had already been primed for immunological memory: dIðsÞ ¼ sλRðsÞ γ þ μ IðsÞ dt dRðsÞ ¼ qðsÞμ þ γ IðsÞ sλRðsÞ μRðsÞ: dt
ð6Þ
System (6) results in basic SIS dynamics, and in the two discrete limits the endemic equilibria become I ¼ ðhsiR0 1Þ=hsiR0 for the homogeneous system, and I ¼ ðhsiR0 1Þ=R0 for the heterogeneous one. These curves are plotted in Fig. 1A,B as dashed lines. Note that sub-model (6) has the same behavior as an extended version of model (2) with a vaccine that has the same mode of action as immunity acquired by natural infection and covers the entire population (Supplementary material, Section S1). This shows that, irrespective of heterogeneity in acquired immunity, the reinfection threshold defines robustly a transmission regime where a vaccine will not be able to eliminate the disease, unless it affords protection that exceeds that conferred by natural exposure (Gomes et al., 2004). These findings can be applied more widely to systems with different routes of transmission as also shown in the Supplementary material, Section S2. 3.3. SIR with distributed acquired immunity and two connectivity groups To illustrate the effects of a form of heterogeneity that affects individuals irrespectively of their infection history, we consider a situation where transmission occurs in a population composed of two risk groups (Rodrigues et al., 2009). Moreover, we assume a positive correlation between propensity to acquire infection and opportunity to transmit to others (Pastor-Satorras and Vespignani, 2001), potentiating superspreading events. This is conceived as a form of contact heterogeneity, modeled by placing each individual in one of two connectivity groups in proportions p1 and p2 , where we conventionally assign subscript 2 to a relatively small group of highly connected (high risk) individuals, while the remaining, with lower contact rates (low risk) are denoted by subscript 1. Heterogeneity in acquired immunity is again represented by the relative susceptibility factor s (0 os o1) with distribution q(s). This subsection refers to Fig. 1C,D.
The fully extended model is given by: dSi ¼ pi μ λ i S i μ S i dt dI i ðsÞ ¼ λi ðqðsÞSi þsRi ðsÞÞ γ þ μ I i ðsÞ dt dRi ðsÞ ¼ γ I i ðsÞ sλi Ri ðsÞ μRi ðsÞ; ð7Þ dt R1 for i ¼ 1; 2, where λi ¼ κ i =hκ i β 0 κ 1 I 1 ðsÞ þ κ 2 I 2 ðsÞds is the force
of infection upon naive individuals of group i, and κ i is the number of contacts made by individuals of group i per unit of time. The basic reproduction number for this model, defined as the spectral radius (dominant eigenvalue) of the next generation operator (Diekmann et al., 1990) and calculated as in Driessche den
(van 2 hκ i β = γ þ μ , where and Watmough, 2002), is given by R ¼ κ = 0
2 2 2 hκ i ¼ p1 κ 1 þ p2 κ 2 and κ ¼ p1 κ 1 þ p2 κ 2 describe the first and second moments of the contact rate distribution in the population. As before, we represent endemic equilibrium infection prevalence in terms of R0, and compare the described trends for qðsÞ following continuous beta distributions and the two discrete limits, all with the same mean hsi. One of the discrete limits corresponds to homogeneous acquired protection and is written as: dSi ¼ pi μ λ i S i μ S i dt dI i ¼ λi ðSi þ σ Ri Þ γ þ μ I i dt dRi ¼ γ I i σλi Ri μRi ; ð8Þ dt for i ¼ 1; 2, where λi ¼ κ i =hκ i β ðκ 1 I 1 þ κ 2 I 2 Þ is the force of infection, and σ is the susceptibility of each recovered subject relative to naive individuals, which is made to coincide with the mean relative susceptibility of recovered subjects in the more general model. The other limit is the all-or-nothing acquired protection written as: dSi ¼ pi μ λ i S i μ S i dt dI i ð0Þ ¼ ð1 σ Þλi Si γ þ μ I i ð0Þ dt dI i ð1Þ ¼ λi ðσ Si þ Ri ð1ÞÞ γ þ μ I i ð1Þ dt dRi ð0Þ ¼ γ I i ð0Þ μRi ð0Þ dt dRi ð1Þ ¼ γ I i ð1Þ λRi ð1Þ μRi ð1Þ; ð9Þ dt where λi ¼ κ i =hκ i βðκ 1 I 1 þ κ 2 I 2 Þ is the force of infection for I i ¼ I i ð0Þ þ I i ð1Þ, and σ is the proportion of individuals who fail to retain protection to subsequent infections. Endemic equilibrium curves are plotted in Fig. 1C,D (in linear and logarithmic scales, respectively). The figure illustrates a setting where 4% of the population has 22.9 times more contacts than the population average (p2 ¼ 0:04 and κ 2 ¼ 22:9), supporting a variance-to-mean ratio of 20. By fixing the value for the mean at hκ i ¼ 1 and usingp1 þ p2 ¼ 1, we get p1 ¼ 0:96 and κ 1 ¼ 0:0875, for the majority of the population with fewer contacts. Throughout the analysis, the effective reproduction number is Re ¼ hsiR0 , consistent with the invariance of the reinfection threshold depicted by the dashed vertical lines. Although characteristics of the contact rate distribution may affect the expression of the reproduction number in terms of various model parameters, the reinfection threshold given by condition Re ¼ 1, remains at R0 ¼ 1=hsi, establishing this critical value as a robust indicator for planning control strategies.
dI i ðsÞ ¼ sλRi ðsÞ γ þ μ I i ðsÞ dt dRi ðsÞ ¼ qðsÞpi μ þ γ I i ðsÞ sλRi ðsÞ μRi ðsÞ; dt
ð10Þ
The endemic equilibria for system (10) are plotted in Fig. 1C,D as dashed lines. This sub-model can also be realized as system (7) extended by a vaccine whose mode of action is identical to immunity acquired by natural infection. As seen for system (2), in system (7) the reinfection threshold defines a region of prevalence in which such a vaccine will not be able to eliminate infection.
variance = 20
0.2 variance = 0 0.1
0.02
variance = 20
0.01 variance = 0 0
proportion infected
3.4. The fragile relationship between R0 and prevalence The comparison between top and bottom panels in Fig. 1 shows that the specification of two risk groups (bottom) critically changes the relationship between endemic prevalence and R0, relative to an implementation where only the mean risk is specified (top). In particular, for the same R0, prevalence is consistently lower under the heterogeneous risk formulation or, conversely, for the same endemic prevalence, R0 is expected to be higher in a heterogeneous population relative to a homogeneous one. Since both the epidemic and reinfection thresholds retain their values on the transmission axis, accounting for realistic heterogeneities implies not only that a population is further above the epidemic threshold (R0 ¼ 1), but also has a greater chance to be above the reinfection threshold (R0 ¼ 2 in the settings of Fig. 1). An equilibrium prevalence of I¼ 0.02, for example, requires R0 ¼1.6 (below the reinfection threshold in Fig. 1A,B) in the homogeneous context and R0 ¼3.3 (above the reinfection threshold in Fig. 1C,D) when two risk groups are considered. This difference is crucially reflected on the estimation of efforts required for control and elimination, as illustrated in Fig. 2B where a universal vaccination programme is implemented for 20 years in the two scenarios (homogeneous in red and heterogeneous in grey). Using parameter values as before while both vaccine and natural infection halve the susceptibility of each individual, we see that the same vaccination programme appears much more effective in the homogeneous setting everything else being equal. Consequently, if realistic heterogeneities in individual risk of infection are not accounted for in the analysis, elimination and control are likely to be unexpectedly challenging. The consequences of neglecting these heterogeneities are less significant when this does not place the study population on the wrong side of critical thresholds. Fig. 2A illustrates a population with prevalence I ¼0.2, which is above the reinfection threshold irrespective of heterogeneity (R0 ¼2.4 by the homogeneous model and R0 ¼101.6 with two risk groups), while Fig. 2C corresponds to a population with prevalence I ¼0.002, positioned below threshold (R0 ¼ 1.1 or R0 ¼1.8, depending on whether the model is homogeneous or heterogeneous). Fig. 3 provides a summary of these results over a range of initial prevalence. According to these analyses, universal coverage vaccination programmes have systematically more impact in populations where individual risk is homogeneous. This deviation can be very large for populations presenting an intermediate initial prevalence of infection, such that R0 is below some critical threshold in the homogeneous scenario and above according to the heterogeneous formulation. For simplicity and maximal effects, this illustration uses a vaccine that perfectly mimics the mode of action of naturally acquired immunity. The phenomenon here depicted, however, that homogeneous models generally
101
0 proportion infected
As done for system (2), we consider a sub-model of system (7), in which the S compartment is bypassed to gain further insight into the reinfection threshold (Gomes et al., 2005):
proportion infected
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
0.002
0.001 variance = 20 variance = 0 0
0
5
10
15
20
time (in years) Fig. 2. Simulation of a universal coverage vaccination programme. Pre-vaccination prevalence is assumed at equilibrium: A, I¼ 0.2; B, I ¼0.02; C, I¼ 0.002. The vaccine mimics natural exposure, halving the susceptibility of each individual. The red curve is the post-vaccination dynamics of the model that reduces individual risks to the mean (6), whereas the grey curve results from the model with two risk groups (10) giving a variance of 20 in individual risk. In both cases the mean individual risk is 1.
underestimate intervention efforts required to meet control targets, is relevant across a broad spectrum of health interventions.
4. Discussion Heterogeneity in individual response to infection and in propensity to transmit to others is not contested, and much research is devoted to understanding their effects. However, this is often missing in population-level modeling of endemic infectious diseases. We use a general transmission model to explain how mean field approximations decrease the estimated R0 for a given prevalence, resulting in underestimation of the effort required for disease control, especially when this approximation inappropriately places transmission estimates below critical thresholds. The theoretical framework constructed here can be adapted to adjust global prevalence datasets for specific endemic infectious diseases. This will enable the identification of invariant transmission thresholds while simultaneously estimating risk distributions and setting-specific transmission parameters (Smith et al., 2005; Gomes et al., 2012) to inform future policies. In the case of tuberculosis, Dye and Ravignione (2013) have recently proposed a broader view of disease risk, moving control from the use of specific vaccines and drugs to a wider context of health promotion and development along with other infectious and non-infectious diseases. Referring to our Fig. 1, this is a change in perspective from intervening on the vertical axis by combating disease directly, to intervening on the horizontal axis by acting on the ensemble of epidemiological and socio-economic factors that affect transmission. Implicitly, this is placing the accurate identification of transmission thresholds at the core of global health policy, making our results very timely and open to applications in many endemic pathogen systems (Klepac et al., 2015).
102
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
prevalence post−vaccination
−1
10
−2
10
−3
10
−3
10
−2
10
−1
10
prevalence pre−vaccination
vaccine effectiveness
100 80 60 40 20 0 −3 10
−2
10
−1
10
prevalence pre−vaccination Fig. 3. Effectiveness of a vaccination programme at 20 years post-vaccination and equilibrium post-vaccination. The vaccine has the same mode of action as immunity induced by natural infection, reducing the force of infection upon naive individuals by 50%. A, Prevalence post-vaccination vs prevalence pre-vaccination. B, Vaccine effectiveness calculated as 1 IðTÞ=Ið0Þ, where T is time since the vaccination programme started and Ið0Þ is the prevalence when the programme starts. Solid lines were calculated 20 years from the beginning of the intervention, while dotted were obtained when the system re-equilibrates after much longer time; grey lines were generated by the heterogeneous model, while red come from the respective mean field approximation (homogeneous). The prevalence of infection at the reinfection threshold, according to each model realization, is marked by a vertical dashed line.
Acknowledgments We thank Professor Antonio Coutinho for constructively challenging the practical usage of earlier versions, urging us to aim for finer concepts. We also thank two anonymous reviewers for constructive comments. MGMG is grateful to Instituto Gulbenkian de Ciência for hosting and to FCT (IF/01346/2014) and CAPES (Science without Borders) for funding. CR is supported by Fundação para a Ciência e Tecnologia, UID/MAT/04561/2013.
Appendix A. Supplementary material Supplementary data associated with this article can be found in the online version at http://dx.doi.org/10.1016/j.jtbi.2016.01.029.
References Anderson, R.M., May, R.M., 1991. Infectious Diseases of Humans: Dynamics and Control. Oxford University Press, United Kingdom. Anderson, R.M., Medley, G.F., May, R.M., Johnson, A.M., 1986. A preliminary study of the transmission dynamics of the human immunodeficiency virus (HIV), the causative agent of AIDS. IMA J. Math. Appl. Med. Biol. 3, 229–263.
Andreasen, V., 2011. The final size of an epidemic and its relation to the basic reproduction number. Bull. Math. Biol. 73, 2305–2321. Ball, F., 1985. Deterministic and stochastic epidemic models with several kinds of susceptibles. Adv. Appl. Probab. 17, 1–22. Becker, N.G., Bahrampour, A., Dietz, K., 1995. Threshold parameters for epidemics in different community settings. Math. Biosci. 129, 189–208. Boylan, R.D., 1991. A note on epidemics in heterogeneous populations. Math. Biosci. 105, 133–137. Clancy, D., Pearce, C.J., 2013. The effect of population heterogeneities upon spread of infection. J. Math. Biol. 67, 963–987. Diekmann, O., Heesterbeek, J.A.P., Metz, J.A.J., 1990. On the definition and the computation of the basic reproduction ratio R0 in models for infectious diseases in heterogeneous populations. J. Math. Biol. 28, 365–382. Dowdy, D.W., Golub, J.E., Chaisson, R.E., Saraceni, V., 2012. Heterogeneity in tuberculosis transmission and the role of geographic hotspots in propagating epidemics. PNAS 109, 9557–9562. Dye, C., Ravignione, M., 2013. Weigh all TB risks. Nature 502, S13. Fox, J.P., Elveback, L., Scott, W., Gatewood, L., Ackerman, E., 1971. Herd immunity: basic concept and relevance to public health immunization practices. Am. J. Epidemiol. 94, 179–189. Gomes, M.G.M., White, L.J., Medley, G.F., 2004. Infection, reinfection, and vaccination under suboptimal immune protection: epidemiological perspectives. J. Theor. Biol. 228, 539–549. Gomes, M.G.M., White, L.J., Medley, G.F., 2005. The reinfection threshold. J. Theor. Biol. 236, 111–113. Gomes, M.G.M., Aguas, R., Lopes, J.S., Nunes, M.C., Rebelo, C., Rodrigues, P., Struchiner, C.J., 2012. How host heterogeneity governs tuberculosis reinfection. Proc. Biol. Sci. 279, 2473–2478. Gomes, M.G.M., Lipsitch, M., Wargo, A.R., Kurath, G., Rebelo, C., et al., 2014. A missing dimension in measures of vaccination impacts. PLoS Pathog. 10 (3), e1003849. Heffernan, J.M., Smith, R.J., Wahl, L.M., 2005. Perspectives on the basic reproductive ratio. J. R. Soc. Interface 2, 281–293. Hyman, J.M., Stanley, E.A., 1988. Using mathematical models to understand the AIDS epidemic. Math. Biosci. 90, 415–473. Jones, J.H., Handcock, M.S., 2003. An assessment of preferential attachment as a mechanism for human sexual network formation. Proc. Biol. Sci. 270, 1123–1128. Katriel, G., 2012. The size of epidemics in populations with heterogeneous susceptibility. J. Math. Biol. 65, 237–262. Klepac, P., Funk, S., Hollingsworth, T.D., Metcalf, C.J.E., Hampson, 2015. Six challenges in the eradication of infectious diseases. Epidemics 10, 97–101. Leigh Brown, A.J., et al., 2011. Transmission network parameters estimated from HIV sequences for a nationwide epidemic. J. Infect. Dis. 204, 1463–1469. Lloyd-Smith, J.O., Schreiber, S.J., Kopp, P.E., Getz, W.M., 2005. Superspreading and the effect of individual variation on disease emergence. Nature 438, 355–359. Mossong, J., Hens, N., Jit, M., Beutels, P., Auranen, K., et al., 2008. Social contacts and mixing patterns relevant to the spread of infectious diseases. PLoS Med. 5 (3), e74. Murphy, B.M., Singer, B.H., Anderson, S., Kirschner, D., 2002. Comparing epidemic tuberculosis in demographically distinct heterogeneous populations. Math. Biosci. 180, 161–185. Novozhilov, A.S., 2012. Epidemiological models with parametric heterogeneity: deterministic theory for closed populations. Math. Model. Nat. Phenom. 7, 147–167. Pastor-Satorras, R., Vespignani, A., 2001. Epidemic dynamics and endemic states in complex networks. Phys. Rev. E 63, 066117. Proscham, F., Sethuraman, J., 1976. Stochastic comparison of order statistics from heterogeneous populations, with applications in reliability. J. Multivar. Anal. 6, 608–616. Riley, R.L., 1974. Airborne infection. Am. J. Med. 57, 466–475. Rodrigues, P., Margheri, A., Rebelo, C., Gomes, M.G.M., 2009. Heterogeneity in infection rates can explain high reinfection rates. J. Theor. Biol. 259, 280–290. Rothman, K.J., 2012. Epidemiology: an Introduction, Second ed. Oxford University Press, United Kingdom. Smith, D.L., Dushoff, J., Snow, R.W., Hay, S.I., 2005. The entomological inoculation rate and its relation to the prevalence of Plasmodium falciparum infection in African children. Nature 438, 492–495. Smith, D.L., Mckenzie, F.E., Snow, R.W., Hay, S.I., 2007. Revisiting the basic reproduction number for malaria and its implications for malaria control. PLoS Biol. 5, e42. van den Driessche, P., Watmough, J., 2002. Reproduction numbers and subthreshold endemic equilibria for compartmental models of disease transmission. Math. Biosci. 180, 29–48. Vaupel, J.W., Yashin, A.I., 1985. Heterogeneity's ruses: some surprising effects of selection on population dynamics. Am. Stat. 39, 176–185. Vynnycky, E., Nagelkerke, N., Borgdorff, M.W., van Soolingen, D., van Embden, J.D.A., Fine, P.E.M., 2001. The effect of age and study duration on the relationship between 'clustering' of DNA fingerprint patterns and the proportion of tuberculosis disease attributable to recent transmission. Epidemiol. Infect. 126, 43–62. Ypma, R.J.F., Altes, H.K., van Soolingen, D., Wallinga, J., van Ballegooijen, W.M., 2013. A sign of superspreading in tuberculosis: highly skewed distribution of genotypic cluster sizes. Epidemiology 24, 395–400.
Contents lists available at ScienceDirect
Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/yjtbi
A theoretical framework to identify invariant thresholds in infectious disease epidemiology M. Gabriela M. Gomes a,b,c,n, Erida Gjini d, Joao S. Lopes d, Caetano Souto-Maior d, Carlota Rebelo e a
CIBIO-InBIO, Centro de Investigação em Biodiversidade e Recursos Genéticos, Universidade do Porto, Portugal Instituto de Matemática e Estatística, Universidade de São Paulo, Brazil c Liverpool School of Tropical Medicine, United Kingdom d Instituto Gulbenkian de Ciência, Oeiras, Portugal e Departamento de Matemática e Centro de Matemática, Aplicações Fundamentais e Investigação Operacional, Universidade de Lisboa, Portugal b
H I G H L I G H T S
Host heterogeneity modifies the relationship between R0 and disease prevalence. Neglecting heterogeneity results in underestimated efforts to meet control targets. Invariant transmission thresholds are robust indicators for control planning.
art ic l e i nf o
a b s t r a c t
Article history: Received 14 April 2015 Received in revised form 17 January 2016 Accepted 20 January 2016 Available online 8 February 2016
Setting global strategies and targets for disease prevention and control often involves mathematical models. Model structure is typically subject to intense scrutiny, such as confrontation with empirical data and alternative formulations, while a less frequently challenged aspect is the widely adopted reduction of parameters to their average values. Focusing on endemic diseases, we use a general transmission model to explain how mean field approximations decrease the estimated R0 from prevalence data, while threshold phenomena – such as the epidemic and reinfection thresholds – remain invariant. This results in an underestimation of the effort required to control disease, which may be particularly severe when the approximation inappropriately places transmission estimates below important thresholds. These concepts are widely applicable across endemic pathogen systems. & 2016 Elsevier Ltd. All rights reserved.
Keywords: Endemic infection Epidemic threshold Reinfection threshold Heterogeneity Global health
1. Introduction Traditional epidemiologic research classifies individuals by similarity of stipulated characteristics and conceives models for disease distribution in a population that is compartmentalized (Rothman, 2012). In infectious diseases, hosts move between compartments according to processes, such as disease progression and transmission, describing dynamic patterns of disease in the population (Anderson and May, 1991). The logic of compartmental study design motivates compartmental models calibrated by average parameters: the “mean field” approximation. Of particular interest from such models are threshold parameters, such as the n Corresponding author at: CIBIO-InBIO, Centro de Investigação em Biodiversidade e Recursos Genéticos, Universidade do Porto, Portugal. E-mail address: [email protected] (M.G.M. Gomes).
http://dx.doi.org/10.1016/j.jtbi.2016.01.029 0022-5193/& 2016 Elsevier Ltd. All rights reserved.
basic reproduction number R0, and the effective reproduction number Re, which have become instrumental in the design and evaluation of control strategies (Diekmann et al., 1990; Becker et al., 1995; Anderson and May, 1991). Threshold parameters, generally dependent on model structure, can only be estimated indirectly and subject to assumptions whose impact is often unrecognized (Heffernan et al., 2005). Here we assess the behavior of transmission threshold parameters under heterogeneity in individual capacity to acquire immunity upon infection and in the rate at which contacts are made with other individuals. Populations are undoubtedly composed of individuals that differ in their propensity to acquire infections and in their potential to transmit to others. Overall, variation results from a mix of genetic and environmental factors, including social and physical aspects. Due to demographic processes, differential selection is likely to occur within groups, resulting in patterns of disease for the population as a whole that differ from the
98
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
expectation (Proscham and Sethuraman, 1976; Ball, 1985; Vaupel and Yashin, 1985). For this reason, variation in individual characteristics and group structure must be considered throughout the entire research process that establishes immunization practices, as recognized almost 50 years ago in a classical paper by Fox et al. (1971). Today there is an extensive literature dealing with mathematical analyses of how herd effects of vaccination programs help to reduce incidence and prevalence of infections among the human population. Until recently, most of this work has assumed that within each compartment, individuals are equally susceptible to infection, equally infectious if infected, and equally active in their social contacts. Although deviations of heterogeneous systems from mean field approximations have been well characterized mathematically for epidemic scenarios (Boylan, 1991; Andreasen, 2011; Clancy and Pearce, 2013; Katriel, 2012; Novozhilov, 2012), translations and applications in epidemiology research, especially with regards to endemic pathogens, remain limited. In this paper, we recall a few exceptions that are worth retaining before proceeding with general considerations on the relationship between threshold parameters and disease prevalence. Implication for disease control and elimination are discussed in Section 4.
studies implicate the droplet nucleus mechanism in the transmission of tuberculosis, measles, influenza and many other infections of the respiratory tract (Riley, 1974). Over the years several authors have independently challenged the mean field approximation in TB transmission systems. Murphy et al. (2002) evoked different genetic susceptibility distributions, demographic factors, and transmission intensity, to explain the wide variation in endemic TB levels between countries. Gomes et al. (2012) analyzed published data relating the reinfection proportion among recrudescences and TB incidence to estimate the decomposition into low and high risk groups that could best reproduce the worldwide dataset. Similarly, an independent study based in the city of Rio de Janeiro, Brazil (Dowdy et al., 2012), estimated relative contributions of high and low risk groups compatible with global estimates. Heterogeneity in propensity of infectious individuals to cause secondary cases has also been postulated to explain the observed skewness in genotypic cluster size distributions of M. tuberculosis (Vynnycky et al., 2001; Ypma et al., 2013) and outbreak dynamics in a number of directly transmitted diseases (Lloyd-Smith et al., 2005).
3. Endemic diseases and threshold parameters 2. Heterogeneity and transmission routes 2.1. Sexual transmission Since the beginning of AIDS epidemic and the identification of the aetiological agent in the 1980s, it has been evident that heterogeneity in individual sexual behaviors must be considered in the interpretation of epidemiological data (Anderson et al., 1986; Hyman and Stanley, 1988). The notion is so consensual that most textbooks on infectious disease epidemiology written since then have a major chapter or section dedicated to heterogeneity in sexually transmitted diseases. Much research has been devoted to estimating such contact networks in diverse settings and with different methods in order to accurately predict HIV transmission dynamics (Leigh Brown et al., 2011; Jones and Handcock, 2003). 2.2. Vector transmission Relatively recently, Smith et al. (2005) developed a mathematical framework that allows the estimation of heterogeneity in malaria infection rates from the relationship between parasite rates and entomological inoculation rates in multiple populations. The model was fitted to data from more than 90 localities, a distribution of individual susceptibilities was estimated, and the implications for disease control were discussed for this vectortransmitted disease. Notably, this heterogeneity was later found to greatly enlarge the range of R0 for malaria transmission (Smith et al., 2007), and entered the discussion of control strategies aiming at R0 o 1. 2.3. Airborne transmission The study of heterogeneity in airborne diseases has been less explored since parameter estimation becomes more difficult, although heterogeneity in social contacts, for example, has started to be recognized (Mossong et al., 2008). A prominent case in the topic is tuberculosis (TB), a disease transmitted by droplet nuclei carrying Mycobacterium tuberculosis secreted from patients with pulmonary TB when coughing, sneezing, spitting and other respiratory acts. These tiny droplets disperse throughout the air of enclosed spaces, such as rooms and buildings, and transmission occurs when they are inhaled by another person. Epidemiologic
The approaches used by Smith et al. (2005) in malaria, and Gomes et al. (2012) in tuberculosis, have a similar strategy of fitting a relatively simple model simultaneously to epidemiological data from multiple populations to infer common underlying heterogeneities. This enables the simultaneous estimation of key threshold parameters and the exploration of implications for disease control in specific settings (Gomes et al., 2014). We proceed by defining two threshold parameters of interest in the context of the mean field susceptible-infected-recovered (SIR) model, and then by analyzing the effects of adding heterogeneity incrementally. 3.1. Mean field SIR To define the threshold parameters and their behavior under individual heterogeneity, we rely on a minimal SIR model, which describes the rates of change in the proportions of the population that are susceptible (S), currently infected (I), and recovered (R) by the system of differential equations: dS ¼ μ λS μ S dt dI ¼ λðS þ σ RÞ γ þ μ I dt dR ¼ γ I λσ R μR; dt
ð1Þ
where λ ¼ βI is the per capita rate of infection (force of infection) given an effective contact rate, β, between susceptible and infected individuals, γ is the rate at which individuals recover from infection, μ is the birth rate, here assumed equal to the death rate, and σ is the relative susceptibility of recovered individuals with respect to those who have never encountered the pathogen. The system has two threshold parameters: 1) The epidemic threshold (Becker et al., 1995), above which epidemics occur and infection remains endemic, is given by R0 ¼ 1, where R0, defined as the average number of secondary infections generated by one infectious individual in a totally susceptible population, is given by R0 ¼ β = γ þ μ ; 2) The reinfection threshold (Gomes et al., 2004, 2005), above which epidemics occur and infection remains endemic in a population where all individuals have previously experienced
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
99
0
0.6
−2
10 0.4 0.2
−4
0
0
10
10
1
10
0
10
proportion infected
0.06
0.04 −2
10 0.02
0
proportion infected
proportion infected
0.8
proportion infected
10
1
−4
0
1
2
3
4
5
0
basic reproduction number, R
0
0.5
1
susceptibility
0
10
10
1
10 basic reproduction number, R
0
Fig. 1. Endemic equilibrium as a function of R0 and contact heterogeneity. Each colored curve assumes a different distribution of acquired immunity as represented on the right, with the same mean susceptibility hsi ¼ 0:5, but different variance. A, B Contact rates are distributed uniformly [model (2)]. C, D Contact rate distribution assumes some individuals have more contacts (higher risk) than others [model (7) with p1 ¼ 0:96, κ 1 ¼ 0:0875 and p2 ¼ 0:04, κ 2 ¼ 22:9, supporting a variance-to-mean ratio of 20]. Vertical lines mark the transmission thresholds: epidemic threshold R0 ¼ 1 (solid); reinfection threshold Re ¼ ⟨s⟩R0 ¼ 1 (dashed). Colored dashed lines represent the cases where all individuals have previously experienced the pathogen and become sensitized. Left and right panels represent the same curves in linear and logarithmic scales, respectively. The remaining parameters were set at: γ¼ 1 yr 1 and μ¼ 1/70 yr 1.
the pathogen and become partially immunized, is Re ¼ 1, where Re is the effective reproduction number given by Re ¼ σ R0 . The endemic infection prevalence at the equilibrium of the mean field SIR model is represented in terms of R0 by the solid red line in Fig. 1A,B (in linear and logarithmic scales, respectively), and the epidemic and reinfection thresholds are represented by the vertical solid and dashed lines, respectively. 3.2. SIR with distributed acquired immunity To assess the behavior of these thresholds under heterogeneity in the level of immunity acquired by exposure to the pathogen, we formulate an extended model in which the relative susceptibility of individuals in the R compartment, denoted by s and defined between 0 and 1, is described by a probability density function, R1 qðsÞ with mean ⟨s⟩ ¼ 0 qðsÞsds. In this sub-section, which refers to Fig. 1A,B and to the Supplementary material, Section S1, we analyse the model extension just described, formally written as: dS ¼ μ λS μS dt dIðsÞ ¼ λðqðsÞS þ sRðsÞÞ γ þ μ IðsÞ dt dRðsÞ ¼ γ IðsÞ λsRðsÞ μRðsÞ; ð2Þ dt R 1 where λ ¼ β 0 IðsÞds is the force of infection upon naive individuals. The basic reproduction number for this model has the same dependence on basic parameters as the mean field approximation, R0 ¼ β = γ þ μ , and so does the epidemic threshold, since the reinfection process, by definition, does not affect R0. As for the effective reproduction number, it is now written as Re ¼ hsiR0 , showing that the expression for the reinfection threshold is also unaffected by heterogeneity in reinfection since only the mean susceptibility of recovered individuals appears in the formula. We represent endemic equilibrium infection prevalence in terms of R0, and compare the described trends assuming that qðsÞ is a beta
distribution with mean hsi ¼ σ . We choose a beta distribution because it fits the criterion of being defined in the interval (0,1) and allows a spectrum of shapes ranging from homogeneous to all-or-nothing (central panels in Fig. 1). In addition, we consider two limiting discrete scenarios: a homogeneous implementation where all previously infected individuals recover with the same factor of partial protection, qðsÞ ¼ δðs σ Þ; for heterogeneous we adopt a polarized (all-or-nothing) implementation where some recovered individuals retain total protection after infection while the remaining retain none, qðsÞ ¼ ð1 σ ÞδðsÞ þ σδðs 1Þ. The results using beta distributions are obtained numerically, whereas the discrete limits are treated analytically as described below. The model with homogeneous acquired protection is written as (1) above. This previously analysed system (Gomes et al., 2004, for example) admits a disease-free equilibrium for any value of the parameters. This state is stable when R0 o 1, and undergoes a transcritical bifurcation at R0 ¼ 1, losing its stability and giving rise to a stable endemic equilibrium when R0 4 1. The endemic equilibrium for system (1) can be solved explicitly, and in particular the proportion infectious is given by: I¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 ffi hsi R0 γ þ μ μ γ þ μ þ hsi2 R0 γ þ μ μ þ γ þ μ 2hsi R0 γ μ þ μ γ þ μ : 2hsiR0 γ þ μ
ð3Þ The model with all-or-nothing acquired protection is written as: dS ¼ μ λS μ S dt dIð0Þ ¼ ð1 σ ÞλS γ þ μ Ið0Þ dt dIð1Þ ¼ λðσ S þ Rð1ÞÞ γ þ μ Ið1Þ dt dRð0Þ ¼ γ Ið0Þ μRð0Þ dt dRð1Þ ¼ γ Ið1Þ λRð1Þ μRð1Þ dt
ð4Þ
where λ ¼ βI is the force of infection for I ¼ Ið0Þ þ Ið1Þ, and σ is the
100
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
proportion of individuals who fail to retain protection to subsequent infections (susceptibility factor 1) while 1 σ is the proportion acquiring total protection (susceptibility factor 0). Solving system (4) for the endemic equilibrium we find that the proportion infectious is given by: I¼
hsiR0 γ þ ðR0 1Þμ
γ þμ þ
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 hsiR0 γ þ ðR0 1Þμ þ γ þ μ 4hsiR0 γ γ þ μ : 2R0 γ þ μ
ð5Þ The solutions given by expressions (3) and (5) are represented in Fig. 1A,B by the solid red and blue lines, respectively, as functions of R0 . Homogeneous protection maintains prevalence levels that are consistently higher than those maintained under the heterogeneous version with the same mean. Two intermediate scenarios, with continuously distributed acquired protection (system (2)), were treated numerically by discretizing the beta distributions and solving the resulting ordinary differential equations in MATLAB (yellow and green in Fig. 1A,B). Noticeably, partial protection upon recovery generates the same reinfection threshold at R0 ¼ 1=hsi in all cases. Following Gomes et al. (2005), further insight into the reinfection threshold is obtained by considering a sub-model of system (2), in which the S compartment is bypassed and all individuals enter the population into the respective R classes as if they had already been primed for immunological memory: dIðsÞ ¼ sλRðsÞ γ þ μ IðsÞ dt dRðsÞ ¼ qðsÞμ þ γ IðsÞ sλRðsÞ μRðsÞ: dt
ð6Þ
System (6) results in basic SIS dynamics, and in the two discrete limits the endemic equilibria become I ¼ ðhsiR0 1Þ=hsiR0 for the homogeneous system, and I ¼ ðhsiR0 1Þ=R0 for the heterogeneous one. These curves are plotted in Fig. 1A,B as dashed lines. Note that sub-model (6) has the same behavior as an extended version of model (2) with a vaccine that has the same mode of action as immunity acquired by natural infection and covers the entire population (Supplementary material, Section S1). This shows that, irrespective of heterogeneity in acquired immunity, the reinfection threshold defines robustly a transmission regime where a vaccine will not be able to eliminate the disease, unless it affords protection that exceeds that conferred by natural exposure (Gomes et al., 2004). These findings can be applied more widely to systems with different routes of transmission as also shown in the Supplementary material, Section S2. 3.3. SIR with distributed acquired immunity and two connectivity groups To illustrate the effects of a form of heterogeneity that affects individuals irrespectively of their infection history, we consider a situation where transmission occurs in a population composed of two risk groups (Rodrigues et al., 2009). Moreover, we assume a positive correlation between propensity to acquire infection and opportunity to transmit to others (Pastor-Satorras and Vespignani, 2001), potentiating superspreading events. This is conceived as a form of contact heterogeneity, modeled by placing each individual in one of two connectivity groups in proportions p1 and p2 , where we conventionally assign subscript 2 to a relatively small group of highly connected (high risk) individuals, while the remaining, with lower contact rates (low risk) are denoted by subscript 1. Heterogeneity in acquired immunity is again represented by the relative susceptibility factor s (0 os o1) with distribution q(s). This subsection refers to Fig. 1C,D.
The fully extended model is given by: dSi ¼ pi μ λ i S i μ S i dt dI i ðsÞ ¼ λi ðqðsÞSi þsRi ðsÞÞ γ þ μ I i ðsÞ dt dRi ðsÞ ¼ γ I i ðsÞ sλi Ri ðsÞ μRi ðsÞ; ð7Þ dt R1 for i ¼ 1; 2, where λi ¼ κ i =hκ i β 0 κ 1 I 1 ðsÞ þ κ 2 I 2 ðsÞds is the force
of infection upon naive individuals of group i, and κ i is the number of contacts made by individuals of group i per unit of time. The basic reproduction number for this model, defined as the spectral radius (dominant eigenvalue) of the next generation operator (Diekmann et al., 1990) and calculated as in Driessche den
(van 2 hκ i β = γ þ μ , where and Watmough, 2002), is given by R ¼ κ = 0
2 2 2 hκ i ¼ p1 κ 1 þ p2 κ 2 and κ ¼ p1 κ 1 þ p2 κ 2 describe the first and second moments of the contact rate distribution in the population. As before, we represent endemic equilibrium infection prevalence in terms of R0, and compare the described trends for qðsÞ following continuous beta distributions and the two discrete limits, all with the same mean hsi. One of the discrete limits corresponds to homogeneous acquired protection and is written as: dSi ¼ pi μ λ i S i μ S i dt dI i ¼ λi ðSi þ σ Ri Þ γ þ μ I i dt dRi ¼ γ I i σλi Ri μRi ; ð8Þ dt for i ¼ 1; 2, where λi ¼ κ i =hκ i β ðκ 1 I 1 þ κ 2 I 2 Þ is the force of infection, and σ is the susceptibility of each recovered subject relative to naive individuals, which is made to coincide with the mean relative susceptibility of recovered subjects in the more general model. The other limit is the all-or-nothing acquired protection written as: dSi ¼ pi μ λ i S i μ S i dt dI i ð0Þ ¼ ð1 σ Þλi Si γ þ μ I i ð0Þ dt dI i ð1Þ ¼ λi ðσ Si þ Ri ð1ÞÞ γ þ μ I i ð1Þ dt dRi ð0Þ ¼ γ I i ð0Þ μRi ð0Þ dt dRi ð1Þ ¼ γ I i ð1Þ λRi ð1Þ μRi ð1Þ; ð9Þ dt where λi ¼ κ i =hκ i βðκ 1 I 1 þ κ 2 I 2 Þ is the force of infection for I i ¼ I i ð0Þ þ I i ð1Þ, and σ is the proportion of individuals who fail to retain protection to subsequent infections. Endemic equilibrium curves are plotted in Fig. 1C,D (in linear and logarithmic scales, respectively). The figure illustrates a setting where 4% of the population has 22.9 times more contacts than the population average (p2 ¼ 0:04 and κ 2 ¼ 22:9), supporting a variance-to-mean ratio of 20. By fixing the value for the mean at hκ i ¼ 1 and usingp1 þ p2 ¼ 1, we get p1 ¼ 0:96 and κ 1 ¼ 0:0875, for the majority of the population with fewer contacts. Throughout the analysis, the effective reproduction number is Re ¼ hsiR0 , consistent with the invariance of the reinfection threshold depicted by the dashed vertical lines. Although characteristics of the contact rate distribution may affect the expression of the reproduction number in terms of various model parameters, the reinfection threshold given by condition Re ¼ 1, remains at R0 ¼ 1=hsi, establishing this critical value as a robust indicator for planning control strategies.
dI i ðsÞ ¼ sλRi ðsÞ γ þ μ I i ðsÞ dt dRi ðsÞ ¼ qðsÞpi μ þ γ I i ðsÞ sλRi ðsÞ μRi ðsÞ; dt
ð10Þ
The endemic equilibria for system (10) are plotted in Fig. 1C,D as dashed lines. This sub-model can also be realized as system (7) extended by a vaccine whose mode of action is identical to immunity acquired by natural infection. As seen for system (2), in system (7) the reinfection threshold defines a region of prevalence in which such a vaccine will not be able to eliminate infection.
variance = 20
0.2 variance = 0 0.1
0.02
variance = 20
0.01 variance = 0 0
proportion infected
3.4. The fragile relationship between R0 and prevalence The comparison between top and bottom panels in Fig. 1 shows that the specification of two risk groups (bottom) critically changes the relationship between endemic prevalence and R0, relative to an implementation where only the mean risk is specified (top). In particular, for the same R0, prevalence is consistently lower under the heterogeneous risk formulation or, conversely, for the same endemic prevalence, R0 is expected to be higher in a heterogeneous population relative to a homogeneous one. Since both the epidemic and reinfection thresholds retain their values on the transmission axis, accounting for realistic heterogeneities implies not only that a population is further above the epidemic threshold (R0 ¼ 1), but also has a greater chance to be above the reinfection threshold (R0 ¼ 2 in the settings of Fig. 1). An equilibrium prevalence of I¼ 0.02, for example, requires R0 ¼1.6 (below the reinfection threshold in Fig. 1A,B) in the homogeneous context and R0 ¼3.3 (above the reinfection threshold in Fig. 1C,D) when two risk groups are considered. This difference is crucially reflected on the estimation of efforts required for control and elimination, as illustrated in Fig. 2B where a universal vaccination programme is implemented for 20 years in the two scenarios (homogeneous in red and heterogeneous in grey). Using parameter values as before while both vaccine and natural infection halve the susceptibility of each individual, we see that the same vaccination programme appears much more effective in the homogeneous setting everything else being equal. Consequently, if realistic heterogeneities in individual risk of infection are not accounted for in the analysis, elimination and control are likely to be unexpectedly challenging. The consequences of neglecting these heterogeneities are less significant when this does not place the study population on the wrong side of critical thresholds. Fig. 2A illustrates a population with prevalence I ¼0.2, which is above the reinfection threshold irrespective of heterogeneity (R0 ¼2.4 by the homogeneous model and R0 ¼101.6 with two risk groups), while Fig. 2C corresponds to a population with prevalence I ¼0.002, positioned below threshold (R0 ¼ 1.1 or R0 ¼1.8, depending on whether the model is homogeneous or heterogeneous). Fig. 3 provides a summary of these results over a range of initial prevalence. According to these analyses, universal coverage vaccination programmes have systematically more impact in populations where individual risk is homogeneous. This deviation can be very large for populations presenting an intermediate initial prevalence of infection, such that R0 is below some critical threshold in the homogeneous scenario and above according to the heterogeneous formulation. For simplicity and maximal effects, this illustration uses a vaccine that perfectly mimics the mode of action of naturally acquired immunity. The phenomenon here depicted, however, that homogeneous models generally
101
0 proportion infected
As done for system (2), we consider a sub-model of system (7), in which the S compartment is bypassed to gain further insight into the reinfection threshold (Gomes et al., 2005):
proportion infected
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
0.002
0.001 variance = 20 variance = 0 0
0
5
10
15
20
time (in years) Fig. 2. Simulation of a universal coverage vaccination programme. Pre-vaccination prevalence is assumed at equilibrium: A, I¼ 0.2; B, I ¼0.02; C, I¼ 0.002. The vaccine mimics natural exposure, halving the susceptibility of each individual. The red curve is the post-vaccination dynamics of the model that reduces individual risks to the mean (6), whereas the grey curve results from the model with two risk groups (10) giving a variance of 20 in individual risk. In both cases the mean individual risk is 1.
underestimate intervention efforts required to meet control targets, is relevant across a broad spectrum of health interventions.
4. Discussion Heterogeneity in individual response to infection and in propensity to transmit to others is not contested, and much research is devoted to understanding their effects. However, this is often missing in population-level modeling of endemic infectious diseases. We use a general transmission model to explain how mean field approximations decrease the estimated R0 for a given prevalence, resulting in underestimation of the effort required for disease control, especially when this approximation inappropriately places transmission estimates below critical thresholds. The theoretical framework constructed here can be adapted to adjust global prevalence datasets for specific endemic infectious diseases. This will enable the identification of invariant transmission thresholds while simultaneously estimating risk distributions and setting-specific transmission parameters (Smith et al., 2005; Gomes et al., 2012) to inform future policies. In the case of tuberculosis, Dye and Ravignione (2013) have recently proposed a broader view of disease risk, moving control from the use of specific vaccines and drugs to a wider context of health promotion and development along with other infectious and non-infectious diseases. Referring to our Fig. 1, this is a change in perspective from intervening on the vertical axis by combating disease directly, to intervening on the horizontal axis by acting on the ensemble of epidemiological and socio-economic factors that affect transmission. Implicitly, this is placing the accurate identification of transmission thresholds at the core of global health policy, making our results very timely and open to applications in many endemic pathogen systems (Klepac et al., 2015).
102
M.G.M. Gomes et al. / Journal of Theoretical Biology 395 (2016) 97–102
prevalence post−vaccination
−1
10
−2
10
−3
10
−3
10
−2
10
−1
10
prevalence pre−vaccination
vaccine effectiveness
100 80 60 40 20 0 −3 10
−2
10
−1
10
prevalence pre−vaccination Fig. 3. Effectiveness of a vaccination programme at 20 years post-vaccination and equilibrium post-vaccination. The vaccine has the same mode of action as immunity induced by natural infection, reducing the force of infection upon naive individuals by 50%. A, Prevalence post-vaccination vs prevalence pre-vaccination. B, Vaccine effectiveness calculated as 1 IðTÞ=Ið0Þ, where T is time since the vaccination programme started and Ið0Þ is the prevalence when the programme starts. Solid lines were calculated 20 years from the beginning of the intervention, while dotted were obtained when the system re-equilibrates after much longer time; grey lines were generated by the heterogeneous model, while red come from the respective mean field approximation (homogeneous). The prevalence of infection at the reinfection threshold, according to each model realization, is marked by a vertical dashed line.
Acknowledgments We thank Professor Antonio Coutinho for constructively challenging the practical usage of earlier versions, urging us to aim for finer concepts. We also thank two anonymous reviewers for constructive comments. MGMG is grateful to Instituto Gulbenkian de Ciência for hosting and to FCT (IF/01346/2014) and CAPES (Science without Borders) for funding. CR is supported by Fundação para a Ciência e Tecnologia, UID/MAT/04561/2013.
Appendix A. Supplementary material Supplementary data associated with this article can be found in the online version at http://dx.doi.org/10.1016/j.jtbi.2016.01.029.
References Anderson, R.M., May, R.M., 1991. Infectious Diseases of Humans: Dynamics and Control. Oxford University Press, United Kingdom. Anderson, R.M., Medley, G.F., May, R.M., Johnson, A.M., 1986. A preliminary study of the transmission dynamics of the human immunodeficiency virus (HIV), the causative agent of AIDS. IMA J. Math. Appl. Med. Biol. 3, 229–263.
Andreasen, V., 2011. The final size of an epidemic and its relation to the basic reproduction number. Bull. Math. Biol. 73, 2305–2321. Ball, F., 1985. Deterministic and stochastic epidemic models with several kinds of susceptibles. Adv. Appl. Probab. 17, 1–22. Becker, N.G., Bahrampour, A., Dietz, K., 1995. Threshold parameters for epidemics in different community settings. Math. Biosci. 129, 189–208. Boylan, R.D., 1991. A note on epidemics in heterogeneous populations. Math. Biosci. 105, 133–137. Clancy, D., Pearce, C.J., 2013. The effect of population heterogeneities upon spread of infection. J. Math. Biol. 67, 963–987. Diekmann, O., Heesterbeek, J.A.P., Metz, J.A.J., 1990. On the definition and the computation of the basic reproduction ratio R0 in models for infectious diseases in heterogeneous populations. J. Math. Biol. 28, 365–382. Dowdy, D.W., Golub, J.E., Chaisson, R.E., Saraceni, V., 2012. Heterogeneity in tuberculosis transmission and the role of geographic hotspots in propagating epidemics. PNAS 109, 9557–9562. Dye, C., Ravignione, M., 2013. Weigh all TB risks. Nature 502, S13. Fox, J.P., Elveback, L., Scott, W., Gatewood, L., Ackerman, E., 1971. Herd immunity: basic concept and relevance to public health immunization practices. Am. J. Epidemiol. 94, 179–189. Gomes, M.G.M., White, L.J., Medley, G.F., 2004. Infection, reinfection, and vaccination under suboptimal immune protection: epidemiological perspectives. J. Theor. Biol. 228, 539–549. Gomes, M.G.M., White, L.J., Medley, G.F., 2005. The reinfection threshold. J. Theor. Biol. 236, 111–113. Gomes, M.G.M., Aguas, R., Lopes, J.S., Nunes, M.C., Rebelo, C., Rodrigues, P., Struchiner, C.J., 2012. How host heterogeneity governs tuberculosis reinfection. Proc. Biol. Sci. 279, 2473–2478. Gomes, M.G.M., Lipsitch, M., Wargo, A.R., Kurath, G., Rebelo, C., et al., 2014. A missing dimension in measures of vaccination impacts. PLoS Pathog. 10 (3), e1003849. Heffernan, J.M., Smith, R.J., Wahl, L.M., 2005. Perspectives on the basic reproductive ratio. J. R. Soc. Interface 2, 281–293. Hyman, J.M., Stanley, E.A., 1988. Using mathematical models to understand the AIDS epidemic. Math. Biosci. 90, 415–473. Jones, J.H., Handcock, M.S., 2003. An assessment of preferential attachment as a mechanism for human sexual network formation. Proc. Biol. Sci. 270, 1123–1128. Katriel, G., 2012. The size of epidemics in populations with heterogeneous susceptibility. J. Math. Biol. 65, 237–262. Klepac, P., Funk, S., Hollingsworth, T.D., Metcalf, C.J.E., Hampson, 2015. Six challenges in the eradication of infectious diseases. Epidemics 10, 97–101. Leigh Brown, A.J., et al., 2011. Transmission network parameters estimated from HIV sequences for a nationwide epidemic. J. Infect. Dis. 204, 1463–1469. Lloyd-Smith, J.O., Schreiber, S.J., Kopp, P.E., Getz, W.M., 2005. Superspreading and the effect of individual variation on disease emergence. Nature 438, 355–359. Mossong, J., Hens, N., Jit, M., Beutels, P., Auranen, K., et al., 2008. Social contacts and mixing patterns relevant to the spread of infectious diseases. PLoS Med. 5 (3), e74. Murphy, B.M., Singer, B.H., Anderson, S., Kirschner, D., 2002. Comparing epidemic tuberculosis in demographically distinct heterogeneous populations. Math. Biosci. 180, 161–185. Novozhilov, A.S., 2012. Epidemiological models with parametric heterogeneity: deterministic theory for closed populations. Math. Model. Nat. Phenom. 7, 147–167. Pastor-Satorras, R., Vespignani, A., 2001. Epidemic dynamics and endemic states in complex networks. Phys. Rev. E 63, 066117. Proscham, F., Sethuraman, J., 1976. Stochastic comparison of order statistics from heterogeneous populations, with applications in reliability. J. Multivar. Anal. 6, 608–616. Riley, R.L., 1974. Airborne infection. Am. J. Med. 57, 466–475. Rodrigues, P., Margheri, A., Rebelo, C., Gomes, M.G.M., 2009. Heterogeneity in infection rates can explain high reinfection rates. J. Theor. Biol. 259, 280–290. Rothman, K.J., 2012. Epidemiology: an Introduction, Second ed. Oxford University Press, United Kingdom. Smith, D.L., Dushoff, J., Snow, R.W., Hay, S.I., 2005. The entomological inoculation rate and its relation to the prevalence of Plasmodium falciparum infection in African children. Nature 438, 492–495. Smith, D.L., Mckenzie, F.E., Snow, R.W., Hay, S.I., 2007. Revisiting the basic reproduction number for malaria and its implications for malaria control. PLoS Biol. 5, e42. van den Driessche, P., Watmough, J., 2002. Reproduction numbers and subthreshold endemic equilibria for compartmental models of disease transmission. Math. Biosci. 180, 29–48. Vaupel, J.W., Yashin, A.I., 1985. Heterogeneity's ruses: some surprising effects of selection on population dynamics. Am. Stat. 39, 176–185. Vynnycky, E., Nagelkerke, N., Borgdorff, M.W., van Soolingen, D., van Embden, J.D.A., Fine, P.E.M., 2001. The effect of age and study duration on the relationship between 'clustering' of DNA fingerprint patterns and the proportion of tuberculosis disease attributable to recent transmission. Epidemiol. Infect. 126, 43–62. Ypma, R.J.F., Altes, H.K., van Soolingen, D., Wallinga, J., van Ballegooijen, W.M., 2013. A sign of superspreading in tuberculosis: highly skewed distribution of genotypic cluster sizes. Epidemiology 24, 395–400.
