The role of interconnectivity in control of an Ebola epidemic

Several West African countries - Liberia, Sierra Leone and Guinea - experienced significant morbidity and mortality during the largest Ebola epidemic to date, from late 2013 through 2015. The extent of the epidemic was fueled by outbreaks in large urban population centers as well as movement of the pathogen between populations. During the epidemic there was no known vaccine or drug, so effective disease control required coordinated efforts that include both standard medical and community practices such as hospitalization, quarantine and safe burials. Due to the high connectivity of the region, control of the epidemic not only depended on internal strategies but also was impacted by neighboring countries. In this paper, we use a deterministic framework to examine the role of movement between two populations in the overall success of practices designed to minimize the extent of Ebola epidemics. We find that it is possible for even small amounts of intermixing between populations to positively impact the control of an epidemic on a more global scale.


A.1 Parameterization
As stated in the main text we choose an R 0 value of 1.85, which is consistent with estimates of the 2013-2015 West Africa Ebola epidemic ?. We assume that the estimate of R 0 is based on transmission in the absence of hospital and community interventions. Therefore, we select values of β H and β F for which R 0 = 1.85 holds without intervention, i.e.b H =b Q = 0 andb U = 1. This is captured by the white line in Figure 1. For all simulations, the various curves correspond to different relative contributions of transmission from funerals and from undetected individuals to the overall R 0 . All other model parameters are provided in Table 1.

A.2 Derivation of R 0 : non-spatial model
Here, we derive the basic reproductive number, R 0 , for the non-spatial model using the next generation matrix method. The relevant classes are E, I U , I H , I Q , and F . Note that although I Q does not contribute to the force of infection, it is the source of some new individuals into the F class. We only include new infections in, F , the "gains" matrix, which is given by: Next, V is the "loss" matrix, and it is given by: R 0 is the spectral radius of F V −1 , and it is given by A full description of the interpretation of R 0 is included in the main text.

A.3 Model with spatial coupling
We extend our model to account for migration between two distinct subpopulations. To do this mechanistically, we explicitly take into account the migrating subpopulations, thus allowing the ability to incorporate populations of different sizes ?. Every individual has a home population (denoted by superscript x or y), and there is migration between populations of all living individuals (susceptible, exposed, infectious and recovered). We assume there is no migration of individuals in the funeral F and effectively removed D class. Individuals leave their home population at rate ρ and return at rate τ . Infectious individuals may be too sick to migrate and have slightly reduced migration parameters, ρ ≤ ρ andτ ≤ τ , although for simplicity we typically assume ρ =ρ and τ =τ . Infection can only occur when both susceptible and infectious individuals are present in the same subpopulation which results in the following force of infection: For each of the classes the first superscript refers to the home population and the second refers to the population the individuals are currently situated in (i.e., I ji indicates infectious individuals who live in population j that are currently in population i). Parameters have subscripts denoting the population they are associated with. Here, transmission rate β will correspond to the population currently situated in. The population N i I is composed of individuals currently situated in population i who are involved in transmission (S, E, I U , I H , F , and R). Here, I Q is excluded from N I because individuals do not contribute to the force of infection.
Using the non-spatial model as a baseline, the dynamics of the populations follow a similar pattern but with the inclusion of migration between populations (the final two terms in each equation below). In the home population, susceptible individuals are born into the population from individuals only in their subpopulation N xx At rate λ x the susceptible individuals in S xx enter the exposed class E xx and become infectious at rate σ: After becoming infectious, individuals enter the undetected class and either remain there were probabilityb U x or enter the hospital (b Qx ) or quarantine (b Qx ): where γ H , γ Q are the rates at which individuals leave their respective infective class. These rates are identical irrespective of the population.
Analogously to the non-spatial model, individuals leaving the infectious classes are assumed to either succumb to infection (with probabilityd j , where j = {H, Q} is the respective infective class) or recover (with probability 1 −d j ). It is assumed that a fraction f j of deceased patients are provided with funerals, during which each patient remains capable of transmission. Individuals are then buried after α days and enter a effectively removed class, D, in which they are no longer capable of transmission. Additionally, the remaining fraction 1 − f j individuals directly enter the D class. We can describe these processes with the following differential equations: where j = H, Q. As we assume, deceased individuals do not migrate, there is only a single home population for F and D.
Upon migrating to the other population, individuals enter temporary populations which follow similar dynamics to the home populations: where x, y = 1, 2, x = y and j = H, Q.
Parameters are derived from the population where the individuals is currently situated and can differ between subpopulations:b j ,d j (different practices for how cases are handled) and α j , f j (different practices in how burials are handled). Parameters characteristic of the disease, such as γ j and σ, do not change between populations. All parameters can be found in Table 1.

A.4 Derivation of R 0 : spatial model
We used the next generation matrix to compute the basic reproductive number, R 0 , for the spatially coupled model.
We denote the fraction of the susceptible population in population i as S i , such that S 1 + S 2 = 1.

Similar to
the non-spatial model, the relevant classes are E i , I i U , I i H , I i Q , and F j , but now there is a specific classes for each subpopulation (i = {xx, xy, yx, yy}, j = {x, y}). Thus, we find our "gains" matrix to be: The gains matrix is completely independent of migration into or out of the population but does include the fraction of the susceptible population in each population. This term does not appear in the "loss" matrix, V , which is given by: R 0 is the spectral radius of F V −1 , but because of the complexity of the analytic form of R 0 , we exclude it here. Importantly, the fraction of susceptible individuals does not impact the value of R 0 if 0 < S i < 1. We use the next generation matrix to compute R 0 numerically in our simulations.
Under simplifying assumptions we can recover the expression for R 0 for the nonspatial model. For example, excluding demographics (µ = 0) and intervention (b H1 = b H2 =b Q1 =b Q2 = 0), assuming equal migration (ρ = τ ) and equal population sizes (S 11 + S 21 = S 22 + S 12 ), and assuming all other parameters equal between the two populations (i.e. β U 1 = β U 2 , β F 1 = β F 2 ), the basic reproduction number becomes: Relaxing some of these assumptions to consider when the populations are not identical but continue to exclude migration (τ = ρ = 0), then the next generation matrix has two non-zero eigenvalues: the larger of which is the value for R 0 .
Indeed, in the absence of migration (τ = ρ = 0) the two populations act independently and all of the analysis of the non-spatial model is recovered for each sub population. In this case, the value for R 0 of the full system is the eigenvalue associated with one population, whichever is largest. Figure 3 of the main text, there is a critical combination of β H and f H such that R 0 = 1, independent of β U and β F . We therefore investigate the relationship observed analytically and determine the numerical value of these critical values. For simplicity, we find this value in the absence of migration.

As observed in
In order to simplify R 0 in Figure 3 we assume that all individuals enter the I H class so no individuals remain in the undetected class and no individuals enter the hospital class (b H = 1,b U =b Q = 0), thus: To determine the critical combination of β H and f H , we consider "extreme" values of β F and β U . First, if we assume that β F = 0, then Equation A4 becomes: Given the parameters we used, this corresponds to a value of β H = 0.0542 which matches the vertical line in Figure 3. On the other extreme, if we assume that β U = 0 then: