Characterizing the dynamics underlying global spread of epidemics

Abstract

Over the past few decades, global metapopulation epidemic simulations built with worldwide air-transportation data have been the main tool for studying how epidemics spread from the origin to other parts of the world (e.g., for pandemic influenza, SARS, and Ebola). However, it remains unclear how disease epidemiology and the air-transportation network structure determine epidemic arrivals for different populations around the globe. Here, we fill this knowledge gap by developing and validating an analytical framework that requires only basic analytics from stochastic processes. We apply this framework retrospectively to the 2009 influenza pandemic and 2014 Ebola epidemic to show that key epidemic parameters could be robustly estimated in real-time from public data on local and global spread at very low computational cost. Our framework not only elucidates the dynamics underlying global spread of epidemics but also advances our capability in nowcasting and forecasting epidemics.

Fig. 1

Validating the framework in the two-population model. a–c The analytical (red dashed lines) and simulated (gray lines) pdf of $T_{ij}^1$, $T_{ij}^5$, and $T_{ij}^{10}$ for an exemplary influenza pandemic, where the mean generation time T_g is 3.5 days and the initial epidemic doubling time t_d is 5 days. The epidemic origin has a population size of 7 million and is seeded with 10 infections at time 0. The mobility rate w_ij is 5 × 10^–6, 5 × 10^–5, and 5 × 10^–4 per day, which span the realistic range for populations with 1–10 million people in the WAN (Supplementary Fig. 1). d–f Quantile–quantile (Q–Q) plots for the analytical and simulated quantiles of $T_{ij}^1$, $T_{ij}^5$, and $T_{ij}^{10}$ across 100 epidemic scenarios randomly generated from the following parameter space using Latin-hypercube sampling: doubling time t_d and generation time T_g both between 3 and 30 days, seed size s_i between 1 and 100. Each epidemic scenario is coupled with a set of network parameters randomly generated with mobility rate w_ij between 10^–6 and 10^–3 and population size N_i between 0.1 and 10 million. Simulated quantiles in each scenario are compiled using 10,000 stochastic realizations. In the Q—Q plots, deviations from the diagonal indicate discrepancies between the analytical and simulated quantiles. Data points are colored in blue if the number of exportations is n or above with probability 1, and yellow otherwise. Insets show the corresponding histograms of percent error in $E\left[T_{ij}^n\right]$

Fig. 2

Validating the framework in the WAN-SPT. a, b Schema of the hub-effect and continuous seeding. In this example, the epidemic arrives at population k after population j has imported three infections from the epidemic origin, i.e., $T_{ij}^3<T_{ik}^1<T_{ij}^4$. In the absence of continuous-seeding adjustment, infection trees spawned by the second and subsequent importations in population j are ignored. c Basic network properties of the WAN-SPT with Hong Kong as the epidemic origin (WAN-SPT-HK). d–f Q–Q plots for the analytical and simulated quantiles of EATs for all 2308 populations in the WAN-SPT-HK across all 100 epidemic scenarios considered in Fig. 1 (i.e., 230,800 Q—Q plots in total). Insets show the corresponding histograms of percent error in expected EAT. d EATs for all 246 populations in D_i,1 before (red) and after (blue) adjusting for the hub-effect. e EATs for all 1828 populations in D_i,2 before (red) and after (blue) adjusting for continuous-seeding and path reduction; hub-effect has been adjusted for the epidemic origin and all populations in D_i,1. f EATs for the remaining 234 populations in D_i,3 and D_i,4 after adjusting for the hub-effect, continuous seeding and path reduction. Supplementary Figures 3–5 provide analogous results for the WAN-SPT with other major hubs as the epidemic origin

Fig. 3

Validating the framework in the WAN. The epidemic origin is Hong Kong as in Fig. 2. a Analogous to Fig. 2d–f, with populations in D_i,c (c = 1, 2, 3, 4) color-coded. Analytical EATs are computed using NPP superposition as described in the main text (see The WAN analysis in Methods for algorithmic details). b Density of the data points in a to show that nearly all the 230,800 Q—Q plots align with the diagonal, which indicates congruence between simulated and analytical EATs. Supplementary Figure 7 provides analogous results with other major hubs in the WAN as the epidemic origin

Fig. 4

Inferring key epidemiologic parameters from surveillance data on global and local spread. Red lines and shades indicate posterior medians and 95% credible intervals of parameter estimates, respectively. a Case study of the 2009 influenza A/H1N1 pandemic in Greater Mexico City. The basic reproductive number R₀ is inferred from the observed EATs for the 12 countries seeded by Mexico as formulated in Balcan et al. Blue circles and error bars indicate the R₀ estimates and their 95% confidence intervals in Balcan et al. assuming that the pandemic started in La Gloria on 11, 18, or 25 February 2009. b Case study of the 2014 Ebola epidemic in Montserrado and Margibi, Liberia. The top panel shows the weekly number of confirmed and probable Ebola cases (bars) and the fitted epidemic curve based on parameters estimated from surveillance data up to 21 September 2014. The middle panel shows retrospective real-time estimates (i.e., nowcasting) of reporting proportion, where the x-axis indicates the date of inference. The bottom panel shows retrospective real-time forecasts of the time to the next international case exportation, with and without adjusting for air travel restrictions started in August 2014. Circles and bars indicate the medians and 99% range of forecasts, respectively. Blue horizontal lines indicate the international case exportations occurred on 20 July and 19 September, 2014. Methods and Supplementary Fig. 10 provide more details and sensitivity analysis