Introduction

Define mobility as the collection of processes through which individuals change their current location. Under such a general definition, human mobility is multi-modal and happens in colossal volumes. We give here a very brief and partial review of mobility, focusing on its link with the spread of infectious pathogens of humans. We present a modelling technique, metapopulations, that can be used when considering such issues. The presentation builds on previous work of the author and collaborators, in particular Arino (2009). Note that we assume here that the reader is familiar with the basic steps involved in the study of a model in mathematical epidemiology.

Duration of mobility

The duration of mobility is the length of time an individual spends away from their place of residence during a mobility event; it varies from minutes to years. A simple arbitrary classification follows. Long-term mobility involves durations longer than six months and typically concerns relocation (immigration, move) or long term displacement from the place of residence (refugees). Six months is chosen as a cut point because public health authorities typically require that individuals affiliate with them if they spend this amount of time in their midst. Medium-term mobility involves events lasting from 1 to 6 months, typically related to domestic and foreign migrant workers, with individuals often maintaining two places of residence. Short-term mobility ranges from over 1 day to 1 month and corresponds to travel away from and back to the place of residence for work or leisure. Finally, instantaneous mobility involves events lasting no more than one day. These events can be daily commute but also work-related travel.

Range of mobility

We call range of mobility the distance away from the place of residence that an individual travels for a mobility event. Ranges vary from metres (local) to thousands of kilometres (intercontinental). Human mobility has evolved greatly over the course of human history. For many years, long range mobility was mainly linked to trade. Most humans were mobile only very locally (within a village or to neighbouring villages). It is only since the beginning of the 20th century that travel for leisure has become more prevalent and that as a consequence, more and more humans undertake short-duration long-range trips.

Coexistence of ranges and durations for a given individual

Different durations and ranges of mobilities can occur simultaneously for a given individual. For instance, a migrant worker can be at the same time undertaking a long-term and long-range mobility event while commuting daily to work from their temporary place of residence.

Time to achieve a range

Another important factor in mobility is the time it takes to achieve mobility on a certain range. For instance, long-range continental and even intercontinental mobility has been going on for centuries. However, as technology evolves, the time it takes to achieve a certain range decreases.

Migrations

Migration usually results in a change of the place of residence. Migration has always existed; it is through migration that earth was populated. The International Organization for Migration (IOM) defines a migrant as a person living in a country that is not their country of birth and under this definition, estimates that there are more than 200 million migrants around the world today, with 60% of these migrants living in developed countries. For instance, from 1980 to 2015, there were 7.4 million new Permanent Residents to Canada, out of an estimated total population in 2015 of 35.8 millions (Government of Canada, 2015). There is also a sizeable migration from rich countries; e.g., the French population registered at a foreign embassy was 1.7 million on December 31st, 2015 (Ministère des Affaires étrangères et du Développement international, 2016). Finally, refugees constitute an important factor of international migration, particularly relevant in the context of public health as they represent populations often in precarious health conditions. UNHCR estimates there were 15 million refugees (mid-2015), 1.2 million asylum-seekers (end 2014), 38 million internally displaced persons (end 2014) and 10 million stateless persons.

Medium-term mobility

In today's globalized world, an important source of mobility is that of migrant workers, who spend extensive periods of time working in a location that is not their primary place of residence. It is estimated that in China alone, there were 273.95 million migrant workers in 2014 (National Bureau of Statistics of the People's Republic of China, 2015). Agricultural workers are another source of temporary migrance, often involving temporary work in richer countries. For instance, there are about 45,000 people per year coming to Canada to work in farms (Government of Canada, 2017).

Short-term and instantaneous mobility

In China, the Chunyun period is a 40 days movement period around Chinese New Year. Estimates are it involved 3.6 billion people-trips in 2014 (Wang, Liu, Mao, Hu, & Gu, 2014). The yearly pilgrimage to Mecca involves an average 2.5 million people (including 1.7 million non-Saudis) (Khan et al., 2010). Instantaneous mobility, such as one's daily commute to work, occurs in numbers that are difficult to evaluate but can safely be estimated to represent tens to hundreds of billions of people-kilometres each day.

Fast spatial spread – the case of SARS

The severe acute respiratory syndrome (SARS) propagated worldwide in 2002–2003. Although it caused only 8273 cases in 28 countries, SARS illustrated the capacity for modern transportation methods to facilitate the rapid and global spread of an infectious pathogen. SARS also inflamed imaginations because its case fatality ratio was high (it caused 775 deaths, i.e., a CFR of 9.4%). The spatial extent and temporal organization of that spread is shown in Fig. 1.

Interactions between health systems – polyomyelitis

Public health systems are typically administered at political geography levels: counties, provinces, states or countries. However, pathogens ignore borders and thus mobility-induced interconnections between administrative units have consequences. One striking example is that of poliomyelitis in Nigeria in 2004. Because of a vaccine scare that started in 2002 (Raufu, 2002), authorities in the northern Nigeria state of Kano interrupted polio vaccination in 2004 (Olusanya, 2004). This occurred in a context where neighbouring countries had polio vaccine coverage that was not ideal; see Fig. 2a.

Vaccinations resumed in July the same year (Fleck, 2004), but between the official interruption and the longer lasting public distrust, the damage was done: cases spiked in Nigeria (Cheng, 2008) but also in close countries such as Cameroon (Endegue-Zanga et al., 2016). The effect is clearly visible in Fig. 2b.

Metapopulation models – general introduction

From the discussion above, formulating mathematical models to help understand the spatio-temporal spread of infectious diseases is very important. Many different methods are available. We focus here on metapopulation models, but stress that understanding complex phenomena such as these requires a multi-faceted approach and the interested reader would benefit greatly from looking into other alternatives. For partial differential equations models, see for instance (Rass & Radcliffe, 2003) and the references therein. Also, the early paper (Bailey, 1980) is worth reading. Concerning the spread of diseases on networks, a standard reference is Newman (2002).

In metapopulation models, space is discrete and consists of individual spatial units called patches. Patches range in spatial extent from a house or neighbourhood all the way up to countries or even continents. There are several types of metapopulation models. We focus here on models with explicit movement, in which each patch contains a dynamical system describing the dynamics of the disease in the populations present within the patch and a coupling function describes the explicit movement of individuals between patches. Other model types, not considered here, assume that there is an implicit coupling between patches and do not incorporate movement explicitly.

Assume geographical locations are elements of a set P with cardinality |P| (the number of patches). Each patch is a vertex in a directed graph (digraph). In each patch, there are individuals from a variety of compartments in the sense of Jacquez and Simon (1993). A compartment could contain for instance individuals susceptible to the disease, another individuals infected with the disease; there could be compartments for the different species affected by the disease, etc. We denote C the set of all compartments that exist in the general system and Cp⊆C those compartments present in patch p∈P.

Individuals in a given compartment may move between patches, with mcqp the rate of movement of individuals in compartment c∈C from patch p∈P to patch q∈P∖{p}. This defines a digraph Gc with arcs Ac for each compartment c∈C, with an arc from p to q if mcqp>0 and no arc otherwise. Because there are |C| compartments in the overall system (whether they are present in each patch or not), each pair (p,q) of patches is joined in each direction by at most |C| arrows, so the general context is that of a multi-digraph.

In turn, each patch p∈P is endowed with a dynamical system describing the evolution of numbers or densities of individuals in the compartments in Cp. Usually, all the patches are equipped with the same type of system; this is what we describe here. Also, we assume here that movement is instantaneous between patches and that individuals do not change compartments as they move, so that, for instance, an individual leaves their patch of origin and arrives in the destination patch while susceptible. Some models have assumed that infection can happen during transport, e.g., (Arino, Sun, & Yang, 2016); others have considered the role of travel time, e.g., Liu, Wu and Zhou (2008). Such models are beyond the scope of this introductory review.

The general form of a metapopulation model is then, ∀c∈C and ∀p∈P,(1)xcp′=fcp(Xp)+∑q∈Pmcpqxcq,where Xp=(x1p,…,x|C|p)T∈ℝ|C| is the distribution of individuals of the different compartments present in the patch and(2)mcpp=−∑q∈P∖{p}mcqp,is used to make the notation more compact. The function fcp(Xp) indicates that processes are localized to the patch where they take place. This excludes for example the model in Arino and van den Driessche (2003), wherein birth occurs in the patch of residence but depends on the total resident population of that patch. If it is important to take into account individuals from the same compartment present in other patches, we define Xc=(xc1,…,xc|P|)T∈ℝ|P| and use fcp(Xp,Xc). However, to simplify the discussion, we assume the form in (1) from now on.

The general form (1) can be used to summarize many different types of metapopulation models, not only ones for the spread of infectious pathogens. (Metapopulations have a rich history in their originating field of Ecology; see, e.g., Hanski and Gilpin (1997).) To focus on epidemiology, we let U⊊C and I⊊C be uninfected and infected compartments, respectively, with U∩I=∅ and U∪I=C. Let skp and iℓp be the number or density of individuals in uninfected and infected compartments, respectively, with k∈U, ℓ∈I and p∈P. Then for all k∈U, ℓ∈I and p∈P, (1) can be written(3a)skp′=fkp(Sp,Ip,Np)+∑q∈Pmkpqskq(3b)iℓp′=gℓp(Sp,Ip,Np)+∑q∈Pmℓpqiℓq,where Sp=(s1p,…,s|U|p), Ip=(i1p,…,i|I|p) and Np=∑u∈Usup+∑ℓ∈Iiℓp.

Behaviour of the total population

Consider the total population Np=Sp+Lp+Ip+Rp in (5). We have(6)Np′=bp−dpNp+∑X∈{S,L,I,R}∑q∈PmXpqXq.

Writing X∈{S,L,I,R}, let b=(b1,…,b|P|)T,N=(N1,…,N|P|)T,X=(X1,…,X|P|)T∈ℝ|P| be vectors and d=diag(d1,…,d|P|), MX be |P|×|P|-matrices, with(7)MX=(−∑q∈PmXq1mX12mX1|P|mX21−∑q∈PmXq2mX2|P|mX|P|1mX|P|2−∑q∈PmXq|P|).

Then (6) can be written as(8)N′=b−dN+∑X∈{S,L,I,R}MXX.

It is useful to derive some properties of the movement matrix (7). For a matrix M, define σ(M)={λ∈ℂ;Av=λv,v≠0}, the spectrum of M, τ(M)=max{ℜ(λ);λ∈σ(M)} the spectral absissa of M and ρ(M)=max{|λ|;λ∈σ(M)} the spectral radius of M. Furthermore, we use the same notation as Berman and Plemmons (1994) for matrices and vectors: M≥0 if mij≥0 for all i,j; M>0 if M≥0 and ∃(i,j) such that mij>0; M≫0 if mij>0 for all i,j.Lemma 1Let

Mc

be a

|P|×|P|-movement matrix for compartment

c∈C. Then the following hold true:1.−Mc

is a singular M-matrix.2.τ(Mc)=0.3.τ(Mc)∈σ(Mc).4.One of the left eigenvectors associated to the eigenvalue

τ(Ms)=0

is the vector

1|P|T=(1,…,1).5.If additionally, Mc

is irreducible, then

τ(Mc)=0

has multiplicity 1, 1|P|T

is (to a multiple) the only strongly positive left eigenvector associated with

Mc

and there is a strongly positive eigenvector

v*

corresponding to

τ(Mc).6.Let

D>0

be a diagonal matrix. If

D≫0, then

D−Mc

is a nonsingular M-matrix and

(D−Mc)−1>0. If

Mc

is irreducible, then

D−Mc

is a nonsingular M-matrix and

(D−Mc)−1≫0.PROOF. Points 1–5 are shown for example in Arino and Portet (2015). To prove 6, first consider the case D≫0. Define d¯=minp∈P{dpp}>0. Then −Mc≤d¯I−Mc≤D−Mc. From Fiedler (2008, Theorem 5.2.5), d¯I−Mc is an M-matrix. Since τ(Mc)=0, using a “spectrum shift” (Horn & Johnson, 2013, Problem 1.2.P8), all eigenvalues of d¯I−Mc have real parts larger than d¯, so d¯−Mc is a nonsingular M-matrix. In turn, Fiedler (2008, Theorem 5.1.1.4) implies that D−Mc is a nonsingular M-matrix and Fiedler (2008, Theorem 5.1.1.11) allows to conclude. Suppose now that Mc is irreducible. Let d¯=maxp∈P{dpp}>0. Then D−Mc is irreducibly diagonally dominant with all columns k∈P such that dkk=d¯ satisfying the strict diagonal dominance requirement. As a consequence, Varga (2004, Theorem 1.11) implies that D−Mc is nonsingular and inverse positivity follows from Fiedler (2008), Theorem 5.2.10).

Note also that for each compartment, −Mc is a nonsymmetric Laplacian matrix. This implies that additional results are available about the spectrum of Mc; see, e.g., Agaev and Chebotarev (2005). We characterize the relations between movement matrices as follows. We say movement is equal for all compartments if MS=ML=MI=MR=:M; movement is similar for all compartments if the zero/nonzero patterns in all matrices are the same but not the entries. Finally, movement is dissimilar across compartments if there are two movement matrices with different zero/nonzero patterns.

Case of movement equal for all compartments. Equation (8) reduces to(9)N′=b−dN+MN,giving the equilibrium solution N*=(d−M)−1b, which, by 6 in Lemma 1, exists if all death rates are positive or M is irreducible. Assume either of these properties holds. To show that N* attracts all solutions of (9), consider Df, the Jacobian of f(N):=b−dN+MN. We have Df=M−d. Since d is nonsingular, we have τ(M−d)<0, so limt→∞N(t)=(d−M)−1b. This is summarized as follows.Theorem 2Suppose that movement is equal for all compartments and that all death rates are positive or the (common) movement matrix

M

is irreducible. ThenN*:=limt→∞N(t)≫0.

It is actually quite “difficult” to be in a situation where Image 1: if b≫0 and d−M is a nonsingular M-matrix, then strong positivity of N* is guaranteed. Thus, except if d=0 or b≥0 with some zero entries, N*≫0.

Case of unequal movement. When movement is not equal for all compartments, the problem is more complicated and, to the best of our knowledge, still open in general. To operate within a manageable context, assume movement rates are similar for all compartments. This situation might arise for instance in the case where infectious individuals are less prone to travel. LetM¯=[minX∈{S,L,I,R}mXpq]pq,p≠qM¯=[maxX∈{S,L,I,R}mXpq]pq,p=qandM¯=[maxX∈{S,L,I,R}mXpq]pq,p≠qM¯=[minX∈{S,L,I,R}mXpq]pq,p=q.

Using these matrices in (8) gives the following inequality:b−dN+M¯N≤N′≤b−dN+M¯N.here, M¯ and M¯ are obtained from movement matrices but are not movement matrices in the sense of Lemma 1; in particular, their column sums are not all zero. It follows that the upper bound grows unbounded and the lower bound does not guarantee nonnegativity.

Note that this does not preclude conducting an analysis of the system. Indeed, some authors have considered models with movement similar for all compartments or even dissimilar across compartments (as might arise if infectious individuals are not allowed some movements – as has been the case with HIV positive individuals travelling to some countries, for instance). However, in this case, results such as Theorem 2 are not readily available and the analysis typically relies on properties of the system embedded in the patches.

The disease free equilibrium (DFE)

Assume (5) is at equilibrium with Lp=Ip=0 for p∈P. Then Φp=0 and for all p∈P,0=bp−dpSp+νpRp+∑q∈PmpqSSq0=−(νp+dp)Rp+∑q∈PmpqRRq,which we must solve for Sp,Rp. Here again, we write the problem in vector form. We have0=b−dS*+νR*+MSS*0=−(ν+d)R*+MRR*,where S*,R*,b∈ℝ|P| and d,ν=diag(ν1,…,ν|P|),MS,MR are |P|×|P|-matrices. We have0=−(ν+d)R*+MRR*⇔(MR−ν−d)R*=0,so the solution R*=0 is unique if MR−ν−d is invertible. This is similar to the proof of Theorem 2 and uses Lemma 1. We thus find that limt→∞R(t)=0 provided that ν+d≫0 or Mr irreducible. Provided this is true, at the DFE, there holds that L*=I*=R*=0 and b−dS+MSS*=0, i.e., S*=(d−MS)−1b. To summarize, provided all death rates are positive or MR and MS irreducible, then the DFE takes the form(S*,L*,I*,R*)=((d−MS)−1b,0,0,0).

Note that S*=N* at the DFE.

Computation of R0 and local stability of the DFE

Use the next generation method (van den Driessche & Watmough, 2002) with infected variables L and I ordered asXI={L1,…,L|P|,I1,…,I|P|}.

Write the dynamics of the subsystem of (5) in the variables XI as XI′=F−V, with F=(Φ1,…,Φ|P|,0,…,0)T and V having all other flows (with a negative sign). Differentiating F with respect to XI, first note that whenever k≠p,∂Φp∂Lk=∂Φp∂Ik=0,since there are no contacts between individuals not located in the same patch. Thus,DF=(diag(∂Φ1∂L1,…,∂Φ|P|∂L|P|)diag(∂Φ1∂I1,…,∂Φ|P|∂I|P|)00).

Suppose that Φp=βpSpIp, i.e., ηp=0. Then ∂Φp/∂Lp=0 and ∂Φp/∂Ip=βpSp. At the DFE, ∂Φp/∂Ip=βpNp*. Suppose that Φp=βpSpIp/Np, i.e., ηp=1. Then ∂Φp/∂Lp=−βpSpIp/Np2 and ∂Φp/∂Ip=βpSp/Np−βpSpIp/Np. Thus, at the DFE, ∂Φp/∂Lp=0 and ∂Φp/∂Ip=βp. Thus,∂Φp∂Lp|DFE=0and∂Φp∂Ip|DFE=βp(N*)1−ηp.

As a consequence,(10)F=DF|DFE=(0diag(∂Φ1∂I1|DFE,…,∂Φ|P|∂I|P||DFE)00)=:(0F1200).

On the other hand,V=DV|DFE=(diagp(εp+dp)−ML0−diagp(εp)diagp(γp+dp)−MI),where the notation diagp(zp)=diag(z1,…,z|P|) is used for simplicity. Thus, using the formula for the inverse of a 2×2 block matrix (Horn & Johnson, 2013), Section 7.3,V−1=((diagp(εp+dp)−ML)−10V˜21−1(diagp(γp+dp)−MI)−1),where(11)V˜21−1=(diagp(εp+dp)−ML)−1diagp(εp)(diagp(γp+dp)−MI)−1.

So we obtain the next generation matrixFV−1=(0F1200)(V˜11−10V˜21−1V˜22−1)=(F12V˜21−1F12V˜22−100),where V˜ij−1 is the block ij in V−1. So R0=ρ(F12V˜21−1). It is difficult to go further without more information. However, given parameter values, R0 is easy to compute numerically. In terms of the local asymptotic stability of the DFE, (van den Driessche & Watmough, 2002), Theorem 2 gives the following result.Theorem 3Define

R0

for the

|P|-SLIRS (5) as(12)R0=ρ(F12V˜21−1),with

F12

and

V˜21−1

defined by (10) and (11), respectively. Then the DFE(S*,L*,I*,R*)=((d−MS)−1b,0,0,0)is locally asymptotically stable if

R0<1

and unstable if

R0>1.

Global stability of the DFE

Several methods are available to consider global stability properties of the DFE in metapopulation models and have been used successfully in many cases. Note that it has also been shown that metapopulation models can exhibit backward bifurcations (Arino, Ducrot, & Zongo, 2012), so proving the global asymptotic stability of the DFE will not always be possible.

Existence and stability of an endemic equilibrium point (EEP)

The problem of existence and stability of an endemic equilibrium point where the disease is present in all patches is more complex and cannot yet be attacked in full generality. Some interesting examples are those in Iggidr et al., 2012, Li and Shuai, 2009.

Inheritance of dynamical properties

Most of the problems mentioned above are solved by exploiting the fact that metapopulations are linearly coupled copies of constituting units whose behaviour is known. This poses a very interesting problem, the inheritance problem. Consider model (3) with no movement,(13a)skp′=fkp(Sp,Ip,Np)(13b)iℓp′=gℓp(Sp,Ip,Np),where k∈U, ℓ∈I and p∈P. Think of (13) as representing any model for a single population found in the literature. Within a given patch p, suppose some dynamical properties are known about (13). Which of these properties also hold for the general system with mobility (3)? It is relatively easy to show that such things as existence and boundedness of solutions is preserved. As shown in this article, the same is true to some extent of existence and local asymptotic stability of the DFE. As we have just seen, global asymptotic stability is not as easy to settle. Techniques such as the one in Li and Shuai (2010) allow to envision a general treatment of this issue when a Lyapunov function is known for the isolated model.

Metapopulation specific behaviour – mixed equilibria

A problem essentially dual to that of inheritance consists in seeking properties that hold for the metapopulation system (3) but not for uncoupled population level model (13). A prototypical example here is the problem of mixed equilibria, first considered in Arino and van den Driessche (2003). Consider a metapopulation such as (3) in which, when there is no movement, i.e., in the case of (13), some patches are at the DFE while others are at an EEP. It is shown in Arino (2009) that for such mixed equilibria to persist in the presence of movement, movement must be dissimilar between uninfected and infected compartments.

Deterministic simulations

It is worth noting that numerical simulations of metapopulation models such as (5) are simple to set up and not extremely onerous on modern computers for systems with a reasonable number of patches. Sample code in R and MatLab is provided as electronic appendices to the present paper to illustrate this.

Stochastic simulations

When dealing with problems related to the fast spatial spread of a pathogen such as SARS (Section 4.1), systems of ordinary differential equations allow to investigate an average behaviour but do not capture the range of potential outcomes that occur in an actual “real life” crisis. For this reason, stochastic simulations of metapopulation models have been carried out; see, e.g., Arino et al., 2011, Colizza and Vespignani, 2008.