Theoretical basis to measure the impact of short-lasting control of an infectious disease on the epidemic peak

Abstract

Background

While many pandemic preparedness plans have promoted disease control effort to lower and delay an epidemic peak, analytical methods for determining the required control effort and making statistical inferences have yet to be sought. As a first step to address this issue, we present a theoretical basis on which to assess the impact of an early intervention on the epidemic peak, employing a simple epidemic model.

Methods

We focus on estimating the impact of an early control effort (e.g. unsuccessful containment), assuming that the transmission rate abruptly increases when control is discontinued. We provide analytical expressions for magnitude and time of the epidemic peak, employing approximate logistic and logarithmic-form solutions for the latter. Empirical influenza data (H1N1-2009) in Japan are analyzed to estimate the effect of the summer holiday period in lowering and delaying the peak in 2009.

Results

Our model estimates that the epidemic peak of the 2009 pandemic was delayed for 21 days due to summer holiday. Decline in peak appears to be a nonlinear function of control-associated reduction in the reproduction number. Peak delay is shown to critically depend on the fraction of initially immune individuals.

Conclusions

The proposed modeling approaches offer methodological avenues to assess empirical data and to objectively estimate required control effort to lower and delay an epidemic peak. Analytical findings support a critical need to conduct population-wide serological survey as a prior requirement for estimating the time of peak.

Background

The influenza A (H1N1-2009) pandemic began in early 2009, and rapidly spread worldwide. Mathematical epidemiologists characterized the epidemic and provided key insights into its dynamics from the earliest stages of the pandemic [1]. The transmission potential was quantified shortly after the declaration of emergence [2–6], while statistical estimation and relevant discussion of epidemiological determinants were underway before substantial numbers of cases were reported in many countries [1].

Prior to the pandemic, many countries issued the original pandemic preparedness plans and guidelines, aiming to instruct the public and to advocate community mitigation. The goals of the mitigation have been threefold; (a) to delay epidemic peak, (b) to reduce peak burden on hospitals and infrastructure (by lowering the height of peak) and (c) to diminish overall morbidity impacts [7]. To assess these aspects under different intervention scenarios, various modeling studies have been conducted (e.g. [8–10]), most notably, by simulating the detailed influenza transmission dynamics.

Although simulations have aided our understanding of expected dynamics in realistic situations and in different scenarios, analytical methods that objectively determine the required control effort and that make statistical inference (e.g. evaluation of empirically observed delay) have yet to be developed. Focus on epidemic peak (relating to mitigation goals (a) and (b) above) has been particularly understudied. Goal (c), on the other hand, is readily formulated in terms of the so-called final epidemic size. The time delay of a major epidemic (such as that resulting from international border control) has been explored using simplistic modeling approaches [11, 12]; however, the height and time of an epidemic peak involve nonlinear dynamics, rendering analytical approaches difficult. Despite the mathematical complexity, goals (a) and (b) can be more readily understood from empirical data during early epidemic phase than can goal (c), because an explicit understanding of goal (c) in the presence of interventions requires knowledge of the full epidemiological dynamics over the entire epidemic period.

In the present study, we present a theoretical basis from which the impact of an early intervention on the height and time of epidemic peak may be assessed. As a special case, we consider a scenario in which intervention is implemented only briefly during the early epidemic phase (e.g. unsuccessful containment). We employ a parsimonious epidemic model with homogeneously mixing assumption, because nonlinear epidemic dynamics involve a number of analytical complexities. As a first step towards understanding epidemiological factors that influence the epidemic peak, leading to the eventual statistical inference of relevant effects, we seek fundamental analytical strategies to evaluate the impact of short-lasting control on epidemic peak using the simplest epidemic model [13]. For our model to become fully applicable and to more closely match empirical data, a number of extensions are required. We discuss ways by which these extensions can be practically realized.

Methods

Study motivation

We first present our study motivation. During the early epidemic phase of the 2009 pandemic, many countries initially enforced strict countermeasures to locally contain the epidemic. Early intervention includes, but is not limited to, quarantine, isolation, contact tracing and school closure. Nevertheless, once it was realized that a major epidemic was unavoidable, regions and countries across the world were compelled to downgrade control policy from containment to mitigation. Although mitigation also involves various countermeasures (and indeed, mitigation originally intends to achieve the above mentioned goals (a)-(c)), one desires to know the effectiveness of the unsuccessful containment effort. Among its many outcomes, the present study focuses on the height and time of epidemic peak.

The applicability of our theoretical arguments is not restricted to the switch of control policy. In many Northern hemisphere countries, the start of the major epidemic of H1N1-2009 (which may or may not have been preceded by early stochastic phase) corresponds to the summer school holiday period. Adults also take vacation over a part of this period. In addition to strategic school closure as an early countermeasure against influenza [14, 15], school holiday is known to suppress the transmission of influenza [16], mainly because transmission tends to be maintained by school-age children [2, 17–19]. Following this trend, a decline in instantaneous reproduction number has been empirically observed during the summer holiday period of the 2009 pandemic [20]. Transmission resumes once a new semester starts. The effectiveness of the summer holiday period in lowering and delaying the epidemic peak is, therefore, a matter of great interest.

Both questions are addressed by considering time-dependent increase in the transmission rate. Let β be the transmission rate per unit time in the absence of an intervention of interest (or during the mitigation phase in the case of our first question). Due to intervention (or school holiday) in the early epidemic phase, β is initially reduced by a factor α (0 ≤ α ≤ 1) until time t₁ (Figure 1A). Though transmission rate abruptly increases at time t₁ when the control policy is eased or when the new school semester starts, we observe a reduced height of, and a time delay in, the epidemic peak compared to the hypothetical situation in which no intervention takes place (Figure 1B). More realistic situations may be envisaged (e.g. a more complex step function or seasonality of transmission), but we restrict ourselves to the simplest scenario in the present study.

Figure 1

A scenario for time-dependent increase in the transmission potential. A. Time dependent increase in the transmission rate. In the absence of intervention (baseline scenario), the transmission rate is assumed to be constant β over time. In the presence of early intervention, the transmission rate is reduced by a factor α (0 ≤ α ≤ 1) over time interval 0 to t₁. We assume that the product αβ ileads to super-critical level (i.e. αR(0) > 1 where R(0) is the reproduction number at time 0), and t₁ occurs before the time at which peak prevalence of infectious individuals in the absence of intervention is observed. B. A comparison between two representative epidemic curves (the number of infectious individuals) in a hypothetical population of 100,000 individuals. R(0) = 1.5, α = 0.90 and t₁ = 50 days. The epidemic peak in the presence of short-lasting control is delayed, and the height of epidemic curve is slightly reduced, relative to the case in which control measures are absent.

Epidemic model

Here we consider the simplest form of Kermack and McKendrick epidemic model [13], formulated in terms of ordinary differential equations. The following assumptions are made: (i) the population is homogeneously mixing, (ii) the epidemic occurs in a population in which the majority of individuals are susceptible, (iii) the time scale of the epidemic is sufficiently shorter than the average life expectancy at birth of the host, and we ignore the background demographic dynamics, (iv) the epidemic occurs in a closed constant population without immigration and emigration again justified based on time scale, and (v) once an infected individual recovers, he/she becomes completely and permanently immune against further infections. Let the numbers of susceptible, infectious and recovered individuals at calendar time t be S(t), I(t) and U(t), respectively. We use the notation U(t) for recovered individuals to avoid confusion with the instantaneous reproduction number at calendar time t, R(t). The population size N remains constant over time (N = S(t) + I(t) + U(t)). The so-called SIR (susceptible-infected-recovered) model is written as

(1)

where R(t) is the instantaneous reproduction number (i.e., the average number of secondary cases generated by a single primary case at calendar time t) and γ is the rate of recovery. Given time-dependent transmission rate β(t) and susceptible population size S(t) at time t, R(t) is assumed to be given by

(2)

Although β(t) will be dealt with as a simple step function in the following analysis, we use the general notation to motivate future analysis of more complex time-dependent dynamics. We assume that an epidemic starts at time 0 with an initial condition (S(0), I(0), U(0)) = (S₀, I₀,U₀) where I₀ = 1 and U₀/N ≈ 0, i.e. an epidemic occurs in a population in which the majority of individuals are susceptible at t = 0. Under this initial condition, we consider two different scenarios for R(t). First, a hypothetical scenario in which no intervention takes place, i.e.

(3)

which is hereafter referred to as the baseline scenario. Second, we consider an observed scenario in which an intervention takes place during the early stage of the epidemic. Let t₁ and t_m,0be calendar times at which the intervention terminates, and at which a peak prevalence of infectious individuals is observed in the absence of intervention, respectively. As mentioned above, we assume that the intervention reduces the reproduction number by a factor α (0 ≤ α ≤ 1) for 0 ≤ t < t₁. For t ≥ t₁, we assume that the transmission rate is recovered to β as in (3).

(4)

We assume t₁ < t_m,0, i.e. we consider a scenario in which transmission rate recovers before the time at which peak prevalence is observed in baseline scenario. We further assume that R(t) > 1 for t < t₁. That is, the efficacy α of an intervention effort (or summer holiday) is by itself not sufficient to contain the epidemic.

To illustrate our modeling approaches, we consider the transmission dynamics of pandemic influenza (H1N1-2009), ignoring the detailed epidemiological characteristics (e.g. pre-existing immunity, realistic distribution of generation time and the presence of asymptomatic infection). The initial reproduction number in the absence of interventions R(0) is assumed to be 1.4 [2]. Given that expected values of empirically estimated serial interval ranged from 1.9 to 3.6 days [2, 5, 21–23], the mean generation time 1/γ is assumed to be 3 days [24, 25].

Our study questions are twofold. First, we aim to quantify the decline in peak prevalence (I(t)/N) due to a short-lasting intervention. The peak prevalence of the intervention scenario is always smaller than that of baseline scenario (see below), and we show that this difference can be analytically expressed. Second, we are interested in the time delay in observing peak prevalence in the presence of intervention. We develop an approximate strategy to quantify the difference in times of peak between baseline and intervention scenarios.

Difference in peak prevalence

We move on to consider estimates of peak prevalence in two scenarios. For mathematical convenience, we use the prevalence of infectious individuals (I(t)/N) to consider the epidemic peak. The peak prevalence of infectious individuals is preceded by peak incidence (γR(t)I(t)/N) by approximately the mean infectious period of 1/γ days. As was realized elsewhere [26], analysis of prevalence is easier than that of incidence. Beginning with two sub-equations of system (1), we have

(5)

Note that R(t) is a function of S(t). Integrating (5) in baseline scenario, we obtain [27]

(6)

A theoretical condition for the observation of peak prevalence at time t_m,0is dI(t_m,0)/dt = 0, or equivalently, R(t_m,0) = 1. As evident from equation (2), this condition satisfies S(t_m,0) = γ/β. The peak prevalence I(t_m,0)/N is then given by [28]

(7)

Note that S₀/R(0)N represents the proportion yet to be infected and S₀ ln R(0)/R(0)N is the proportion removed at time t_m,0. Equation (7) indicates that the peak prevalence of SIR model is determined by the initial condition and the transmission potential R(0). It should be noted that S₀/R(0) can be replaced by γ/β, and thus, I(t_m,0) is independent of initial condition for U₀ = 0 (a special case).

In the intervention scenario, equation (6) with replacement of β by αβ applies for t < t₁.

(8)

which provides another initial condition at time t = t₁ for t ≥ t₁. That is, we can also employ (6) to compute peak prevalence for t ≥ t₁ with initial condition (S(t₁), I(t₁),U(t₁)). Again, a condition to observe peak prevalence at time t_m,1is R(t_m,1) = 1, which gives S(t_m,1) = γ/β. The peak prevalence I(t_m,1)/N of the intervention scenario is given by

(9)

Note that R(0) in the above equation refers to βS₀/γ (i.e. we use R(0) in our baseline scenario to permit an explicit comparison between the two scenarios). Inserting right-hand side of (8) into (9), we obtain

(10)

Consequently, relative reduction in peak prevalence due to intervention within time t₁ is ∈ _α = (I(t_m,0) - I(t_m,1))/N, which can be parameterized as

(11)

Equation (11) indicates that the difference of peak prevalence between the two scenarios is determined by four different factors; the relative reduction in reproduction number α due to the intervention, initial condition at time 0, transmission potential R(0), and fraction of susceptible individuals at time t₁ under the intervention. If the initial condition, the transmission dynamics in the absence of interventions (i.e. R(0), β and γ) and t₁ are known, an estimate of α gives S(t₁), yielding an estimate of ∈ _α .

Delay in epidemic peak

The time to observe peak prevalence is analytically more challenging than the height of peak prevalence, because even an approximate estimate requires an analytical solution to the model (1). We propose a parsimonious approximation strategy which leads to more convenient solutions than those discussed in past studies (e.g. [29]). Substituting I(t) in the first sub-equation of (1) by (1/γ)(dU(t)/dt), we have

(12)

For the baseline scenario (i.e. β(t) = β), integrating (12) from time 0 to t,

(13)

Because U(0)/N ≈ 0,

(14)

Subsequently, the third sub-equation of (1) is rewritten as

(15)

Here we impose another key approximation. Because the quantity βU(t)/γ (≈ R(0)U(t)/N) for influenza (e.g. R(0) = 1.4) tends to be smaller than 1 (especially, before observing epidemic peak), we use a Taylor series expansion, i.e.,

(16)

If the quadratic approximation is inadequate for large R(0), a higher order Taylor polynomial function can be used. Inserting the quadratic approximation into (15), and imposing a further approximation (i.e. S₀≈ N), we obtain

(17)

which appears to be a logistic equation. We use this logistic-form solution instead of the more commonly employed hyperbolic-form solution [29, 30], to illustrate a simpler approximate solution and to demonstrate the problem underlying both solutions. Later, we use a more formal solution (of logarithmic-form) in the intervention scenario, which is numerically identical to the classical hyperbolic-form solution (see below). Assuming that U(0) = U₀> 0, the analytical solution of (17) is

(18)

The derivative of (18) is dU(t)/dt = γI(t). It follows that

(19)

Further differentiation of (19) gives dI(t)/dt, and letting dI(t)/dt = 0, the time to observe peak prevalence is analytically derived. For the logistic equation, the corresponding time has been referred to as the inflection point of the cumulative curve in equation (18) [31]. The inflection point t_m,0to observe peak prevalence is

(20)

which depends on initial condition and transmission characteristics. In the intervention scenario (in which intervention is short-lasting), an identical approach can be taken for t < t₁, replacing β by αβ (or by replacing R(0) by αR(0)). Subsequently, the epidemic peak occurs at t_m,1(> t₁). We take a similar approach to that used in (15) with a computed initial condition (S(t₁), I(t₁), U(t₁)) using (18) and (19). For t ≥ t₁,

(21)

Now we apply an approximation

(22)

It should be noted that, in the above approximation, we include the term exp(βU(t₁)/γ), because U(t) - U(t₁) better satisfies the Taylor series approximation than expanding U(t) alone. Let constants A, B and C be

(23)

Given these constants, we consider

(24)

where z = t - t₁ and V (z) = U(z + t₁) for t ≥ t₁. The initial condition V (0) is V₀ = U(t₁). Writing (24) in integral form, we have [30]

(25)

Past studies have typically assumed hyperbolic-form functions for the analytical solution of (25) [29, 30]. However, we express the solution in logarithmic-form [31, 32], because logarithmic functions are compatible with spreadsheet programs. The logarithmic-form solution reads

(26)

where

(27)

Also,

(28)

Differentiating (26) with respect to z and taking dV (z)/dz = 0, we find the inflection point z_m,1to be

(29)

Replacing z by t, we obtain

(30)

as an approximate solution of the epidemic peak in intervention scenario. The time delay of this peak, imposed by the intervention in the early epidemic phase, τ _α is subsequently calculated as

(31)

using (20) and (30) for the right-hand side. The delay depends on initial condition U₀, the length of intervention t₁ (both of which are apparent from (20) and (30)) and on the efficacy of intervention α (since this quantity influences the initial condition U(t₁) in (30)).

Application and illustration

Empirical analysis of influenza A (H1N1-2009)

Here, we apply the above described theory to empirical influenza A (H1N1-2009) data. Figure 2 shows the estimated number of influenza cases based on national sentinel surveillance in Japan from week 31 (week ending 2 August) 2009 to week 13 (week ending 28 March) 2010. The estimates follow an extrapolation of the notified number of cases from a total of 4800 randomly sampled sentinel hospitals to the actual total number of medical facilities in Japan. The cases represent patients who sought medical attendance and who have met the following criteria, (a) acute course of illness (sudden onset), (b) fever greater than 38.0°C, (c) cough, sputum or breathlessness (symptoms of upper respiratory tract infection) and (d) general fatigue, or who were strongly suspected of the disease undertaking laboratory diagnosis (e.g. rapid diagnostic testing). Although the estimates of sentinel surveillance data involve various epidemiological biases and errors, we ignore these issues in the present study. Prior to week 31, the number of cases was small and the dynamics in the early stochastic phase have been examined elsewhere [17]. We arbitrarily assume that the major epidemic starts in week 31.

Figure 2

Estimated weekly incidence of influenza cases in Japan from 2009-10. The estimates are based on nationwide sentinel surveillance, covering the period from week 31 in 2009 to week 13 in 2010. The estimate follows an extrapolation of the notified number of cases from a total of 4800 randomly sampled sentinel hospitals to the total number of medical facilities in Japan. The case refers to influenza-like illness cases with medical attendance, possibly involving other diseases, but with influenza A (H1N1-2009) dominant among the isolated influenza viruses during the period of interest. Period A corresponds to summer school holiday, followed by autumn semester (period B). Period C covers winter holiday and period D corresponds to winter semester.

It is interesting to observe that the period A in Figure 2 corresponds to that of summer school holiday. Due to reporting delay of approximately 1 week [17], we assume that weeks 31 to 36 inclusive (the latter of which ends on 6 September) reflect the transmission dynamics during the summer school holiday. Subsequently, school opens in September with an epidemic peak in late November (period B), followed by abrupt decline during the winter holiday (period C) and start of winter semester (period D). Among these periods, we focus on the impact of summer holiday (period A), relative to period B, in lowering epidemic peak and delaying the time to observe the peak. More specifically, we estimate the reproduction number R(0) and its reduction α from the data set encompassing weeks 31 to 42. To permit an explicit estimation, we assume that linear approximation holds, as was similarly assumed elsewhere [15]. We assume that the reproduction number is reduced by a factor α from week 31 to 36 due to summer holiday, while the reproduction number recovers to R(0) from week 37 to 42.

Let r₀ be the exponential growth rate of cases per day in the absence of summer holiday. Because our SIR model approximates the generation time by an exponential distribution with mean 1/γ days, the estimator of R(0) is [33, 34]

(32)

Throughout the summer holiday, we assume that the reproduction number is reduced to R _A = αR(0). That is, the growth rate during the summer holiday, r₁, is defined by

(33)

Let the weekly incidence be J _k in week k. During summer holiday period, the conditionally expected value of J_k+1given J _k is

(34)

where Δt is the length of reporting (i.e. 7 days). In week 37, the conditional expectation is

(35)

because, with an initial incidence i _k in week k, we have

(36)

and

(37)

See [35] for more details regarding the derivation of (35), which has been applied to influenza A (H1N1-2009) in other settings [36]. From week 38 to 42, E(J_k+1; J _k ) simplifies to

(38)

Although adding the information of test negative individuals could potentially yield a less biased estimate of r₀[37], we do not have access to this data and so we disregard this issue for now. Assuming that the observed counts of cases are Poisson distributed within each period, the likelihood function to estimate r₀ and α is

(39)

The maximum likelihood estimates of r₀ and α are obtained by minimizing the negative logarithm of (39), and the 95% confidence intervals (CI) are computed by profile likelihood. From the maximum likelihood estimates, we compute the differences in peak prevalence and times to observe the peak between baseline and second scenarios using the SIR model (1). For simplicity, we adopt (S₀, I₀, U₀) = (99998, 1, 1) for the numerical computation and t₁ is assumed to be 50 days (roughly corresponding to the length of period A plus 8 days).

Although our model (1) adopts exponentially distributed generation time, we can partially address the uncertainty of r₀ and α in the parametric assumption of the generation time. That is, we adopt constant generation time (i.e. delta function) as an alternative assumption, which is known to yield a theoretical maximum reproduction number, given identical r₀ and mean generation time [34]. Given r₀, the estimator of R(0) with constant generation time of 1/γ days reads

(40)

As mentioned above, the reproduction number during summer holiday reduces to R _A = αR(0). Accordingly, the growth rate during the summer holiday, r₁ is written as

(41)

Using the above mentioned likelihood (39) and replacing r₁ of exponential assumption by (41), we estimate r₀ and α. It should be noted that the difference between (33) and (41) indicates that the estimates of both r₀ and α depend on the realistic distribution of the generation time [25].

Sensitivity analysis

In addition to the analysis of empirical data, we also examine sensitivity of the height and time of peak prevalence to different values of α and t₁ by numerical simulation. As mentioned above, influenza is our case study, and the default value of R(0) of baseline scenario is 1.4, but we also consider R(0) of 1.2 and 1.6. These ranges are adopted, additionally, because we impose an approximation (16). When obtaining the numerical solutions, initial condition is fixed at (S₀, I₀, U₀) = (99998, 1, 1) for clarity. U₀ = 1 is adopted to prevent U₀ = 0 in (18) and in later equations, and also to select a positive integer value closest to 0 such that U₀/N ≈ 0.

Results

Influenza A (H1N1-2009)

Figure 3A compares observed and predicted numbers of influenza cases in Japan from week 31 to 42. Grey bars represent conditionally expected values during summer holiday, and white bars represent the expected values during autumn semester. The estimated growth rate in the absence of summer holiday, r₀ is 0.048 (95% CI: 0.029, 0.066) per day. Thus, the estimated reproduction number R(0) is 1.14 (95% CI: 1.09, 1.20) which is likely an underestimate (see below).

Figure 3

The impact of summer holiday on the transmission dynamics of influenza A (H1N1-2009). A. Comparison of the observed and predicted weekly counts of the estimated number of influenza cases in Japan from week 31 to 42. Grey bars are conditionally expected values during the summer school holiday, while white bars are the conditionally expected values during autumn semester. B. Numerical solutions of the number of infectious individuals in a hypothetical population of 100,000 individuals with initial condition (S₀, I₀, U₀) = (99998, 1, 1). Baseline scenario is compared against intervention scenario under a short-lasting intervention. For both scenarios, assumed R(0) and mean generation time are 1.14 and 3 days, respectively. In the second scenario, the transmission rate is reduced by a factor α = 0.948 due to summer school holiday from time 0 to t₁ (=50 days). Although heights of peak prevalence do not greatly differ from each other, the time to observe the peak is clearly delayed in the second scenario.

The estimated relative proportion of the reproduction number under mitigation conditions α is 0.948 (95% CI: 0.842, 1.053). Although the confidence limits of α include 1, our argument adopts linear dynamics for longer than 10 weeks (note that the largest number of notifications is seen in week 48), and thus, the reproduction number is conservatively estimated (i.e. over the time period that we examine, the transmission dynamics may become nonlinear); therefore, R(0) is likely underestimated due to the linear approximation adopted in our quantitative illustration). Thus, we believe it is appropriate to regard the reduction in the reproduction number during summer holiday as marginally significant. The estimated reproduction number during school holiday is 1.08.

Even adopting constant generation time of 3 days, r₀ is of a similar order, i.e. 0.048 (95% CI: 0.029, 0.066) per day. The reproduction number R(0), however, is slightly greater (1.15 (95% CI: 1.09, 1.22)) due to its estimator, exp(r₀/γ). α is estimated at 0.942 (95% CI: 0.832, 1.061).

Figure 3B illustrates the number of infectious individuals in a hypothetical population with 100,000 individuals using the estimated α and R(0). In the absence of summer holiday, the epidemic peak would have been observed at Day 171 with I(t_m,0) = 822 cases. In the presence of summer holiday from time 0 to t₁, the peak is delayed to Day 192 with I(t_m,1) = 820 cases. Thus, the estimated α and R(0) do not significantly alter the height of peak prevalence when the effects of summer holiday are included (a difference between the two scenarios of only 2 cases), but the delay effect between the scenarios is as long as 21 days. We do not use our approximate formula for the estimation of time delay in this empirical case study (for reasons explained below). Because R(0) = 1.14 in the absence of intervention, a major epidemic can occur when the initial condition U₀/N is smaller than 1 - 1/α R(0) = 7.4% of the population. Assuming a fixed I₀ = 1, and varying S₀ and U₀ from N - 1 to 0.926N and from 0 to 0.074N, respectively, only slight variations in the reduction of peak prevalence (not greater than 1 case) are observed, but the time delay varies greatly; under a scenario with αR(0) ≈ 1 or with S₀ = 0.926N and U₀ = 0.074N, a possible maximum delay of t₁ = 50 days can be readily envisaged.

Differential peak prevalence

Figure 4A examines the sensitivity of relative peak prevalence to α (i.e. reduction in R(0)) for assumed R(0) of 1.2, 1.4 and 1.6. Because αR(0) < 1 prevents major epidemic during the early epidemic phase, possible ranges of α satisfying αR(0) ≥ 1 vary with R(0). It is worth noting that the relative reduction in peak prevalence is a nonlinear function of α. Largest reduction occurs when α lies within the range 0.90 to 0.95, rather than when α is minimum. Figure 4B examines the sensitivity of relative reduction in the peak prevalence as a function of the time length of intervention t₁ (e.g. the time required to switch control policy from containment to mitigation). Again, to satisfy t₁ < t_m,0, possible ranges of t₁ vary with R(0). The interpretation of Figure 4B is more straightforward than that of Figure 4A. Essentially, the longer the time period of intervention, the larger the potential reduction of prevalence.

Figure 4

Sensitivity analysis of peak prevalence to the efficacy and length of intervention. A. Sensitivity of relative peak prevalence of infectious individuals to reduction in the reproduction number, α. The vertical axis represents I(t_m,1)/I(t_m,0), where I(t_m,0) and I(t_m,1) represent the peak prevalence in the absence and presence of intervention, respectively. The time length of intervention, t₁ is fixed at Day 50. B. Sensitivity of relative peak prevalence of infectious individuals to the time length of the intervention, t₁. α is fixed at 0.90. For both panels, the lines are truncated to satisfy αR(0) ≥ 1 and t₁ <t_m,0where t_m,0is the time to observe peak prevalence in the absence of intervention.

The nonlinear relationship observed in Figure 4A can be explored by combining (11) with approximate solutions. First, because S(t₁) ≈ S(0) exp (-α R(0)U(t₁)/S₀) in the presence of intervention, equation (11) can be expressed in the form

(42)

Second, using an approximate solution of U(t₁) based on logistic equation (18),

(43)

which is a nonlinear function of α. Although the calculation is not shown here due to its mathematical complexity and lack of practical interpretation, the derivative of (43) with respect to α reveals an optimal α yielding the longest delay in Figure 4A.

Delay in epidemic peak

Figure 5A compares epidemic curves of infectious individuals in the absence of intervention between explicit numerical and approximate solutions (i.e. solutions to (1) and (19), respectively). The height of epidemic peak is approximated well for smaller R(0), reflecting the fact that Taylor series expansion is a good approximation to the exponential function. The relationship between R(0) and approximation of epidemic peak height is also analytically expressed. Inserting (20) into (19), the approximate peak prevalence is

Figure 5

Assessing approximation of time to observe epidemic peak. A. Comparison of epidemic curves in the absence of intervention between explicit numerical and approximate solutions. The population size is assumed to be 100,000 with initial condition (S₀, I₀, U₀) = (99998, 1, 1). The generation time is exponentially distributed with mean 1/γ = 3 days. B. Comparison of estimated delay in epidemic peak (imposed by a short-lasting intervention) between explicit numerical and approximate solutions. The assumed length of intervention t₁ is 50 days.

(44)

It should be noted that S₀/R(0) can be replaced by γ/β, and thus, the peak prevalence in the approximated logistic form is independent of the initial condition (indeed, the derivative of a logistic equation is known not to depend on initial condition but rather on the carrying capacity [31]). On the other hand, an explicit solution can depend on initial condition for U(0) > 0, i.e.,

(45)

These two quantities are identical for R(0) = 1, N ≈ S₀ and U₀/N ≈ 0. Otherwise, numerical observation shows that I(t_m,0) > I _approx (t_m,0) for R(0) > 1 and U₀/N ≈ 0. Thus, the smaller the R(0), the better the approximate height of epidemic curve.

While interpretation of the height of peak prevalence is overall straightforward, the time to observe peak prevalence is better captured for larger R(0) (Figure 5A). Clearly, the logistic equation applied to R(0) = 1.2 yields a considerably biased (delayed) time to observation of epidemic peak (with bias longer than 20 days), and thus, we did not apply our approximate solution to the above mentioned case study of influenza A (H1N1-2009). It must be noted that the relationship between R(0) and approximation of epidemic peak in Figure 5A is regulated not only by R(0) but also by the initial condition U₀ in (20); that is, the good agreement of the time of peak between two solutions for R(0) = 1.6 is not only due to R(0) but also to U₀ = 1 in our simulation setting. For suitable values of U₀, good approximations to the time of peak are obtained even for R(0) = 1.2 or smaller (results not shown).

Figure 5B compares the estimated time delay in epidemic peak, induced by an early intervention, between explicit numerical and approximate solutions. For all three R(0) that we investigate, approximation methods result in underestimation of the delay (31). For R(0) = 1.6, the approximate estimate of delay is crudely realized, and its sensitivity to α is close to that of the explicit numerical solution. The approximation is worst for R(0) = 1.2, for which a negative result was yielded for large α. These findings are in concordance with Figure 5A. It should be remembered that approximation of time of epidemic peak can vary with initial condition U₀, indicating that the estimation requires knowledge of U₀ in addition to R(0) and α (also, as we have seen, good approximation depends on judicious choice of R(0) and U₀).

Discussion

In the present study, we have presented fundamental ideas to assess the impact of a short-lasting intervention of an infectious disease on the epidemic peak. As a first step towards explicit evaluation of control effort in lowering and delaying the epidemic peak, we comprehensively described analytical expressions for the difference in the height of, and the time delay in, the epidemic peak gained by intervention, employing a parsimonious homogenous mixing epidemic model. We restricted our focus to a simple step function (Figure 1A) which adequately illustrated the role of summer holiday in lowering and delaying epidemic peak during the influenza (H1N1-2009) pandemic. Our methods show that both the height and the time of epidemic peak can be readily controlled by varying initial conditions at a given point of time at which transmission rate abruptly changes. The proposed method can be extended in future to encompass more realistic multiple steps for the transmission rate. Analytical solution of the simplest form of Kermack and McKendrick model has been undertaken multiple times, and is documented in many key references [27, 29–32, 38, 39] including the original study in 1927 [13]. However, to our knowledge, the present study is the first to offer a theoretical basis on which to assess the impact of an early countermeasure on the epidemic peak employing logistic and logarithmic-form solutions, with a goal to making statistical inferences in the future. In addition, our analytical expressions not only promote epidemiological understanding of epidemic peak traits, but can be worked backwards to determine required control effort to achieve desired goals of height and delay of the epidemic peak. Although manual adjustment of the efficacy of intervention is not practically feasible, our demonstration of a nonlinear relationship between α and decline in the height of epidemic peak should be considered a key issue for optimal intervention management during an early epidemic phase.

Although public health guidelines have tended to advocate control policy switches from containment to mitigation at some point in time, the likely impact of unsuccessful containment to epidemic dynamics has been seldom discussed. Indeed, common illustration of mitigation (embodied in the three distinct goals (a)-(c) described in the Background section) do not account for the time-dependent transmission rate (rather, a guideline illustrates only differential peaks with various reproduction numbers R(0) [7]). Statistical modeling studies tend to focus only on time-dependent decreases in transmission potential, with a focus on the instantaneous reproduction number R(t) < 1 to demonstrate successful control of an infectious disease (e.g. [40–42]). Motivated by these problems, we considered the impact of upward change in the transmission rate during the course of an epidemic. In our case study of influenza H1N1-2009, it was shown that summer vacation did not appreciably lower the height of an epidemic peak, but imposed a substantial time delay. Although our approximate estimation of the delay in epidemic peak was biased towards certain choices of R(0) and the initial number of immune individuals U₀, this does not imply failure of the model. Rather, by means of approximate analytical computation, we have shown in (20) and (30) that the time of epidemic peak critically depends on the fraction of initially immune individuals prior to an epidemic. Although a preceding study instead emphasized the dependence of the time of peak on initial number of infectious individuals [32], other studies with hyperbolic-form solutions emphasize dependency of the time of peak on U₀[27, 31]. We believe that U₀ is more easily quantified than the initial number of infectious individuals I₀ in practical settings. The critical importance of the initial number of immune individuals U₀ is especially highlighted in the influenza A (H1N1-2009) pandemic because of pre-existing immunity [43–46]. Our analytical undertakings indicate that statistical inference of the time delay gained by early control effort will greatly benefit from population-wide seroepidemiological survey. Depending on the quality of approximation for a given combination of R(0) and U₀, one can then decide whether the estimation of delay should be based on analytical or numerical solution.

Despite our motivation to eventually offer a method to estimate the impact of an early intervention, it should be noted that the present study does not account for uncertainty (e.g. confidence interval). Our arguments are based solely on deterministic models, whereas an explicit derivation of the confidence interval requires use of a stochastic Markov jump process [47]. Moreover, potential model extensions are numerous. Relevant factors include, but are not limited to, mobility of host, spatial dynamics, class-age structure (e.g. infection-age dependency), chronological age-structure, social contact patterns, seasonality, strain specificity and immunological dynamics. All of these features would increase the realism, but would greatly complicate analytical inspection, of the model. To illustrate the way forward, we discuss the simplest example of a class-age structured model in the Appendix.

Despite many future tasks to be completed, and our realization that the epidemic peak is vulnerable to heterogeneous patterns of transmission, relevant statistical assessment (including the estimation of R(0), generation time, and incubation period) always starts with a homogeneous modeling assumption [11, 21, 23, 33, 34, 48, 49]. This is particularly true during the early epidemic phase of a pandemic. In this sense, we believe that the present study has successfully offered a methodological avenue to statistically assess empirical data and to assess required control effort to lower and delay epidemic peak.

Conclusions

This study has presented a theoretical basis on which to assess the impact of short-lasting intervention on the epidemic peak of an infectious disease. Employing a homogeneously mixing epidemic model, we derived analytical expressions for the decline in the height of epidemic peak and for the time delay of the peak. Empirical influenza A (H1N1-2009) data were analyzed using a simplistic but practically accessible model, which estimated that the epidemic peak was delayed for 21 days by the summer holiday period in 2009. Approximate logarithmic form solution of the time of epidemic peak appeared to critically depend on initial condition of immune individuals, supporting a need to conduct population-wide serological survey. Despite obvious needs to address various types of heterogeneity, our framework offers a successful methodological avenue to assess relevant empirical data and to advocate required control effort to lower and delay epidemic peak.

Appendix: A way forward

In realistic situations, there is a time delay for newly infected individuals to acquire infectiousness, known as the latent period. This delay is captured by employing the so-called SEIR (susceptible-exposed-infected-recovered) model. In addition to S(t), I(t) and U(t), we consider infected but non-infectious individuals E(t). Assuming that the mean latent period is 1/δ days, the model is written as

(46)

In considering the height of epidemic peak of this system, an equation similar to (6) is derived by employing a Lyapunov function, W(S, E, I) [50]. A different derivation is given elsewhere [51].

(47)

The Lyapunov function is known to yield a constant solution (for any t), thus, dW/dt = 0. This is confirmed by

(48)

At an epidemic peak at time t _m (where we have dE/dt = dI/dt = 0), two obvious conditions apply,

(49)

(50)

Using Lyapunov function in (47), we have

(51)

Writing both sides of (51) in terms of (47), and taking logarithm of both sides, we obtain

(52)

Rearranging (52) leads us to

(53)

which is very close to (6). By varying initial conditions at each time point when a change in transmission rate occurs, equation (53) permits us to measure the impact of public health intervention on the epidemic peak, under the SEIR assumption. Regarding the time of epidemic peak, a straightforward extension of (15) applies, although the analytical solution may be rather complex or may not exist. We have

(54)

S(t) can be replaced as for (14). Expressing E(t) as a function of U(t) requires strong mathematical supports. The problem may be addressed in different ways, but we first consider an analytical solution of dE/dt, i.e.,

(55)

where $\dot{S}(\sigma )$ in the right-hand side is equivalent to dS/dσ for which the derivative of approximation (15) can be used. We have yet to derive a simple analytical solution of the time of epidemic peak from (54), but the above discussion demonstrates that it is at least possible to compute the height of epidemic peak with more E compartments (as was discussed in [52]), using our proposed Lyapunov approach. Nevertheless, it should be remembered that the presence of infection-age dependency in the infectiousness profile is likely to complicate the computation of epidemic peak [53]. The incorporation of other realistic features, however, is greatly aided by the Lyapunov approach. In further studies, we will address epidemic peak with seasonality and age-dependent heterogeneity and in a multi-strain system.