key: cord-0003437-yxrokwgu authors: Granell, Clara; Mucha, Peter J. title: Epidemic spreading in localized environments with recurrent mobility patterns date: 2018-05-03 journal: nan DOI: 10.1103/physreve.97.052302 sha: 05de1af59e01b6b80d9e4c88316caecc9ad4eb56 doc_id: 3437 cord_uid: yxrokwgu The spreading of epidemics is very much determined by the structure of the contact network, which may be impacted by the mobility dynamics of the individuals themselves. In confined scenarios where a small, closed population spends most of its time in localized environments and has easily identifiable mobility patterns—such as workplaces, university campuses, or schools—it is of critical importance to identify the factors controlling the rate of disease spread. Here, we present a discrete-time, metapopulation-based model to describe the transmission of susceptible-infected-susceptible-like diseases that take place in confined scenarios where the mobilities of the individuals are not random but, rather, follow clear recurrent travel patterns. This model allows analytical determination of the onset of epidemics, as well as the ability to discern which contact structures are most suited to prevent the infection to spread. It thereby determines whether common prevention mechanisms, as isolation, are worth implementing in such a scenario and their expected impact. The spreading of infectious diseases is strongly dependent on the networked structure of interactions in the population [1, 2] and on the mobility of individuals [3] [4] [5] [6] [7] . A particularly interesting scenario is that where the structure of the social contacts of the individuals is not completely well mixed nor completely structured, but offers an intermediate level of description. These kinds of models are usually referred to as metapopulation models [8] [9] [10] [11] and were first used in the field of population ecology [12] [13] [14] [15] . In such settings, the nodes of the network represent a population, which is occupied by individuals, and the links of the network represent the migration of individuals from one population to another. This scenario is particularly useful in the study of the spreading of epidemics, given that many real-life patterns of interactions happen in structured, localized populations connected by some degree of migration. The populations usually describe small, local environments (e.g., a city, a college dormitory, or an office, depending on the application) where it is plausible to assume that every individual in the population is able to contact any other individual inside the same location with some probability. The underlying network structure (i.e., the links between subpopulations) describes the mobility patterns of individuals among locations, and can be weighted and/or directed. The problem of modeling such scenarios relies on finding the appropriate level of abstraction to grasp the main macroscopic features of the epidemic spreading process for individuals across the particular environment. The analysis of these oversimplified model abstractions is of utmost importance to separate the effect of single parameters on the incidence of the spreading process, yet allowing an analytical approach that could be used for prediction purposes and to test prevention actions. Traditionally, models for epidemic spreading in metapopulations [8] rely on reaction-diffusion equations to account for the epidemic and mobility dynamics, and assume that (I) individuals diffuse like random walkers through the network and (II) subpopulations with the same number of connections are treated as statistically equivalent [16] , thus smoothing over the actual contact network between individuals. While this approach has been usefully applied in many scenarios [17] , its simplified assumptions do not capture some important real-world features. For instance, analysis of human mobility data reveals that human dynamics are often dominated by recurrent patterns where individuals have memory of the location they come from [18] and are highly likely to return to their original location after a short exploration of the network [3, 19] . The typical exploration of the network mostly consists in visiting frequently a limited number of locations, predominantly performing commutes between home and work locations [20] . Additionally, the traditional assumption of statistical equivalence of subpopulations of the same degree, while allowing for an analytic solution of the invasion threshold, makes it impossible to quantify the outreach of an epidemic in a particular subpopulation of the network. In this work we present a discrete-time Markov-chain model [21, 22] for epidemic spreading in structured populations with a recurrent pattern of migrations between the locations in a bipartite network. The aim of this model is to quantify the extent of a susceptible-infected-susceptible (SIS) epidemic in the scenario where each individual spends most of their time between two locations: their residence (e.g., home or college dormitories) and common destinations where mixing with individuals coming from other residence subpopulations happen (e.g., work places, classes, or other common event spaces). Our goal is to discern which parameters modeling such scenarios control the phase transition of the spreading of a disease. In so doing, we can determine whether typical mechanisms of isolation-such as reducing the mobility of infected individuals-are able to contain the spreading of diseases. That is, whether or not such interventions change the critical properties of the spreading process. The paper is organized as follows. In the next section we introduce the formulation for our model for epidemic spreading in localized environments with recurrent, bipartite travel connections. In Sec. III we show the derivation of the epidemic threshold. In Sec. IV we introduce an isolation mechanism for the infected individuals, and present the consequent formulation. Section V is devoted to the results of our analysis, and Sec. VI offers a discussion that concludes our work. Our metapopulation network model considers two types of locations: residences and common sites. Each residence i has an associated population of n i agents. Individuals associated to a given residence are assumed to interact with one other in a well-mixed fashion. A common location, on the other hand, does not have a fixed population associated to it, thus providing a meeting site for mixing individuals from different residences. The distribution of individuals in common areas is determined by the weighted flows W , with elements w ij defining the probability of an individual associated to residence i to travel to common location j . The flows W define a bipartite network structure: no direct connections between different residences nor between different common areas are considered. The dynamics of the model follow a discrete-time reactiondiffusion process. Every day (for each time step), individuals diffuse through the flows determined by W according to the mobility probability p, causing n i p agents to travel to a common location and n i (1 − p) individuals to remain in their residence sites, for each residence i. Once the individuals are in their new location, they react with the other individuals in the subpopulation (what we call the daytime infection step), meaning a susceptible individual gets infected upon contact with another infected individual with probability β. Then, agents return to their residences and another reaction is performed (nighttime infection step). After, individuals who were infected at the beginning of the time step may recover spontaneously with probability μ. It is important to stress two particularities of this model. First, the daytime infection step takes place both in the common locations and in the residence sites, therefore affecting individuals who did not migrate as well as those who did. Second, what we consider a full time step comprises two infection steps (day and night) and one recovery step. We are interested in calculating the fraction of infected individuals assigned to any residence i for each time step t, ρ i (t), whose time evolution is described by the following equation: The interpretation of Eq. (1) is that the fraction of infected individuals assigned to residence site i at time t + 1 is calculated as the fraction of individuals that were already infected in the previous time step and did not recover, plus those individuals who were susceptible and got infected at the end of the time step according to probability i (t), which is defined as where p is the mobility probability and C is the number of subpopulations defined as common areas. The four terms in Eq. (2) refer, in order, to the fraction of individuals that did not travel and got infected in their residence site in the daytime step; the fraction of people that did travel and got infected in the common site of destination; the fraction of individuals that did travel, did not get infected in the common area of destination but got infected in their residence at the nighttime step; and finally, the fraction of people that did not travel, did not get infected in their residence during the daytime step but got infected in the residence in the nighttime step. The expressions for the probabilities of getting infected in residence site i during daytime, in residence site i during nighttime and in common area j are, respectively, where n i is the size of residence i, W is the bipartite connectivity matrix and W k = C j W kj . D refers to the number of residential sites. The number of individuals that remain in subpopulation i is and the number of individuals moving from residence k to common location j is From Eq. (1) we can calculate the solution of the system in the steady state, by assuming that ρ i (t + 1) = ρ i (t) = ρ i . Under the assumption that near the critical onset of the epidemics the fraction of infected individuals is negligible, we can substitute ρ i = i 1. Eq. (1) then reads Substituting i by its expression in Eq. (2), we obtain Substituting D ☼ i , D i , and C j by their respective expressions in Eqs. (3), (4), (5), we have Applying the approximations ( n i and removing the O( 2 i ) terms, the previous expression reduces to We can express the previous equation in the form of an eigenvector problem, where our new expression is and thus we obtain the classical expression in epidemic spreading [23] : where the entries of the matrix M are Each entry M ik accounts for the average number of contacts between one individual of residence i and all the individuals associated to any residence k during a full day. Indeed, the first term of the right-hand side (r.h.s.) of Eq. (14), accounts for the total average number of contacts among individuals of the same residence, while the second term accounts for the number of interactions that take place at the common locations. Additionally, to investigate the effects of realistic isolation in our setup, we prescribe the mobility probability of infected individuals to be p p, thus effectively reducing their mobility through the network. The parameter that controls the relation between the two mobility rates is what we call the isolation factor γ , being p = γp with 0 γ 1. This prescription changes the formulation introduced in Sec. II as follows. First, the calculation of the number of individuals remaining in their residence (n i→i ) and the number of individuals going from residence k to common location j (n k→j ) need to be adjusted to take into account the two mobility probabilities. Now the probability that an individual remains in its original residential patch is ( Second, the terms D ☼ i (t) and C j (t) use, in the original formulation, ρ i as a proxy of the probability of infection in subpopulation i. This is no longer appropriate when the isolation factor is active, given that the individuals that remain in residence i will no longer be an arbitrary mixing of infected and susceptible individuals. Instead, residence i in the daytime step will mostly be populated by infected individuals as p grows smaller. The correct approach is to calculate the conditional probability for an individual from population i to be in the infected state (I ) given that the individual remains in the population during the daytime (R), which is P (I |R) = P (R|I )P (I ) P (R|I )P (I ) + P (R|S)P (S) Using the new prescription, Eq. (3) reads now Following the same rationale, we rewrite the expression for Eq. (5): (19) Once these changes are introduced, we can calculate the epidemic threshold of the model with isolation following the same procedure we explained in Sec. III. After linearizing our equation and solving the eigenvector problem, we obtain the same expected expression of Eq. (13) , but now the entries of matrix M are Note that when the isolation mechanism is not active (γ = 1), p = p and the previous expression reduces to Eq. (14) . From the previous expression we see that the parameters that are able to shift the onset of the epidemics are the connectivity matrix W , the vector of sizes of the residential subpopulations n, the mobility probability p and the isolation factor γ . In the next section we explore the effects of those parameters in the final output of the epidemic process. To validate our model, we crosscheck the results obtained in the numerical solutions of our analytic model with extensive Monte Carlo simulations. A comparison is depicted in Fig. 1 , where we plot the fraction of infected individuals in the whole system in the steady state ρ as a function of the infectivity (13) . For this plot, the number of subpopulations of type residential is D = 25 and there are C = 5 common sites, with equal-sized residential sites of 100 individuals each. The isolation mechanism is inactive (γ = 1, p = p ), the recovery probability μ = 0.1, and the connectivity matrix is an unweighted, fully connected bipartite network. Epidemic threshold β c as a function of the mobility probability p, for different configurations of the number of residential and common sites. We observe in all of them a non-monotonic behavior of β c which increases up to an optimum value of the mobility parameter (p * ) that makes the epidemic threshold maximum (see dashed lines). Here the recovery probability μ = 0.2, all residential subpopulations are of the same size n = 25, the isolation mechanism is disabled (γ = 1) and the underlying topology is an unweighted fully connected bipartite network. parameter β, for four values of the mobility probability p and with isolation inactive. The correspondence of our analytical results with the Monte Carlo simulations is remarkable for values of the infectivity parameter even beyond the epidemic threshold. To analyze the effect that the mobility probability p has on the epidemic threshold, we plot in Fig. 2 the curves of the critical onset of the epidemic, for different configurations on the number of residential and common sites. Here, we want to highlight an interesting feature: the curve of β c does not have a monotonic behavior, instead there is an optimum value of the mobility probability (p * ), which makes the epidemic threshold maximum. Indeed, we observe that p * will be smaller than 0.5 if the number of residential subpopulations exceeds the number of common sites (D > C); greater than 0.5 in the opposite case (D < C), and exactly 0.5 if the number of residential and common sites are equal (D = C), for the case of a fully connected unweighted topology and for residential sites being of the same size. This happens because p * is the value of the mobility probability that causes all subpopulations in the network to be of the same (or most similar) effective size during the daytime infection step. The physics rationale of this effect can be understood as follows: the critical threshold of the epidemics is dominated by the critical threshold of the largest subpopulation, so the maximum epidemic threshold (minimum spreading of the disease) will be achieved when all populations (both residential and common sites) are of similar size. Note that the same phenomenology has been reported in [22] for mono-partite metapopulation networks. Up to now we have supposed homogeneity in the sizes of the subpopulations of type residential, meaning all entries of vector n are equal. Now we explore what is the effect that heterogeneity will have in the epidemic threshold. To do so, Here D = 10, C = 5 and the total number of agents is N = 100. We vary how individuals are distributed among residential sites, ranging from the homogeneous case (all residential sites are of equal size) to the most heterogeneous setup (all individuals reside in the same node). We see that as the heterogeneity increases, the epidemic threshold gets smaller, and the effect of the optimum mobility p * is diluted. we keep the total number of individuals in the population constant, but we redistribute individuals in such a way that the variance of the size distribution increases monotonically. The results are displayed in Fig. 3 , where we observe that, as heterogeneity increases (higher variance values), values of p * are shifted right and the maximum values of β c are less peaked. This reflects that the uneven distribution of sizes of residential subpopulations affects the critical threshold in such a way that the more heterogeneity the more easy for the epidemic to become endemic. Finally, we analyze the role of the isolation factor on the critical properties of the model. In Fig. 4 we show the epidemic threshold curve as a function of mobility probability p, for different settings of the isolation factor γ . We observe an interesting effect: if the isolation is inactive (γ = 1), we see the increase of the epidemic threshold before p * and the subsequent decrease as reported in Fig. 2 ; but as we decrease γ from 1 to 0 (thus gradually restricting the mobility of the infected individuals) the critical behavior of the system becomes more favorable to the epidemic extinction. As the mobility of the infected individuals is more restricted, the epidemic threshold increases with increasing mobility. In the particular example reported in Fig. 4 , the critical behavior of the epidemic threshold is monotonically increasing for values of γ 0.3. We also observe that the change in the curvature of all the epidemic threshold functions coincides at exactly the expected value of p * , which in this case is 0.5 given that the number of residential and common sites is the same (D = C). Following the behavior observed in Fig. 2 , we also tried the configurations D < C and D > C and obtained that the crossing point of all curves is p * > 0.5 and p * < 0.5 respectively, as expected, for the case of unweighted fully connected bipartite connectivity matrices. For this plot we have used a fully connected unweighted bipartite network consisting of ten residential sites and ten common locations. All residential patches are of size 100 and the recovery probability is μ = 0.1. We observe that the curves cross at exactly the expected value p * = 0.5 given that there are exactly the same number of residences and common sites (see main text for a broader explanation). Summarizing, in this work we have proposed an analytical model to explore the spreading of epidemics in localized environments with non-random, recurrent mobility patterns. The critical properties of the epidemic process have been determined and corroborated by simulations. The results show that the main effect of the recurrent mobility is that the epidemic threshold depends on the mobility probability in a non-monotonic way, presenting an optimal value for which the epidemic is most contained. We also show that restricting the mobility of the infected individuals is an effective mechanism to to delay the critical threshold, specially for high values of the mobility. Importantly, the presented approach allows the appropriate modeling of epidemics on realistic scenarios that include recurrent mobility among bipartite structures, such as university campuses, home-to-work commutes or the spreading of disease in cities. The current formulation of this model is applicable to particular cases which may require locations of heterogeneous sizes, weighted connectivity and different topological structures, and allows determining whether isolation strategies are worth implementing in such specific scenarios. The presented model not only offers analytical insights to the very important problem of epidemic spreading in localized environments but could also become a powerful tool to use in data analysis and policy making. Modeling Infectious Diseases in Humans and Animals Proc. Natl. Acad. Sci. USA Human Mobility and Spatial Disease Dynamics Spatial Ecology: The Role of Space in Population Dynamics and Interspecific Interactions Environmental Intelligence Unit Genetics and Evolution of Metapopulations Proceedings of the 22nd International Symposium on Reliable Distributed Systems