key: cord-0275496-4805jz44
authors: Reis, Elohim Fonseca dos; Masuda, Naoki
title: Metapopulation models imply non-Poissonian statistics of interevent times
date: 2021-06-18
journal: nan
DOI: 10.1103/physrevresearch.4.013050
sha: f233cb29a5e4be3f9f0e6f5b95dcd9e75e63a46f
doc_id: 275496
cord_uid: 4805jz44

Interevent times in temporal contact data from humans and animals typically obey heavy-tailed distributions, and this property impacts contagion and other dynamical processes on networks. We theoretically show that distributions of interevent times heavier-tailed than exponential distributions are a consequence of the most basic metapopulation model used in epidemiology and ecology, in which individuals move from a patch to another according to the simple random walk. Our results hold true irrespectively of the network structure and also for more realistic mobility rules such as high-order random walks and the recurrent mobility patterns used for modeling human dynamics.

Contact networks among individuals crucially impact contagion in populations and networks, including spreading of infectious diseases, information, opinions, and rumors. In addition to the network structure, temporal aspects of contacts and edges in networks, collectively referred to as temporal networks, play a central role in such dynamics [1] [2] [3] [4] . The interevent time (IET), defined as the time between consecutive events, is one such aspect. Event-driven stochastic dynamics models on networks most typically assume Poisson processes for generating events, which yield exponential IET distributions. However, IETs from various human activity data obey non-Poissonian statistics, in particular, heavy-tailed distributions [1, 3, 5, 6] . In most (but not all) cases, heavytailed IET distributions slow down contagion and diffusion in epidemic processes [7] [8] [9] [10] [11] [12] [13] [14] , opinion dynamics [15] [16] [17] [18] , evolutionary game dynamics [19] , cascade processes [20] [21] [22] [23] , and random walks [24] [25] [26] [27] . Several mechanisms can generate heavy-tailed IET distributions, including priority queue models [5, 6, [28] [29] [30] [31] [32] , self-exciting processes [33] [34] [35] , mixture of exponentials [36] [37] [38] , and dynamics of nodal states, where mutual agreement of two nodes produces contact events at a high rate [39] .

However, these models do not explain the genesis of heavy-tailed IET distributions in the presence of mobility of individuals. Human and animal individuals move around to meet others, transmitting information and disease, and then depart from each other. In fact, most of face-to-face contact or proximity data have been collected from mobile individuals, and such data show heavy-tailed IET distributions [25, [39] [40] [41] [42] [43] [44] [45] [46] [47] . Previous modeling studies showed that encountering events of random walkers moving on a two-dimensional plane with heterogeneous attractiveness values generate heavy-tailed IET distributions [48] [49] [50] [51] [52] . These studies crucially assumed that each walker is endowed with a randomly assigned at- * naokimas@buffalo.edu tractiveness value and that walkers tend to slow down when they approach attractive walkers nearby. In the present study, we are interested in analytically accounting for heavy-tailed IET distributions under mobility in a simpler manner, i.e., without introducing the concept of attractiveness or any heterogeneity among walkers.

A prevalent approach to simultaneous modeling of mobility and agent-to-agent interaction is metapopulation models [53] , originally introduced in mathematical epidemiology [54] [55] [56] [57] [58] and ecology [59] [60] [61] . In a metapopulation network, nodes represent subpopulations, such as, households, cities, urban areas, or ecological habitats, and edges represent migration routes. Metapopulation models have been successful in describing, e.g., epidemic processes [2, 13, [62] [63] [64] including the COVID-19 pandemic [65] .

In the present study, we provide accounts of heavytailed IET distributions when individuals move around, meet, and separate according to a standard metapopulation network model, in which the individuals perform simple random walks, and its more realistic variants. We analytically show that the IET distribution for a pair of individuals is a mixture of exponential distributions, which robustly holds true for different network structures and mobility rules. Crucially, mixtures of exponential distributions produce IET distributions that are substantially closer to heavy-tailed distributions observed in empirical data than exponential distributions do.

We consider a metapopulation network with N subpopulations (i.e., nodes). The network may be directed and/or weighted. The individuals that populate the network are assumed to be simple random walkers in continuous time. Different walkers move independently of each other and may produce contact events (events for short) between them only when they are copresent in the same subpopulation.

Specifically, the time t for which each walker waits until it leaves the current subpopulation is independently arXiv:2106.10348v2 [physics.soc-ph] 7 Feb 2022 distributed according to an exponential distribution with rate D. The mobility and contact dynamics are illustrated in Fig. 1 . When the walker leaves the ith subpopulation, it moves to a neighboring jth subpopulation with probability A ij / N =1 A i , where A ij is the weight of edge (i, j). Events between a pair of walkers in the same subpopulation occur according to a Poisson process with rate λ, which epidemic process models on metapopulation networks usually assume. In other words, the distribution of time t to the next event between a pair of walkers in the same subpopulation is given by φ e (t) = λe −λt . When the walkers are in different subpopulations, no event occurs between them.

Consider two walkers visiting the same subpopulation. We derive the probability density of an IET, τ , between an arbitrary pair of individuals, which are random walkers. When the two walkers are copresent in a subpopulation, events (including the case of 0 event) occur between them before either walker leaves the subpopulation. After that, a series of n copresences occurs before the next event occurs (see Fig.1 ).

We let p(τ, n) be the probability of generating an IET of length τ during the nth copresence after either walker leaves the current subpopulation. The probability density of the time at which either walker leaves the subpopulation where the two walkers are copresent is given by φ d (t) = 2De −2Dt . We obtain

By Laplace transforming p(τ, 0), we obtain

.

We denote the distribution of time t between consecutive copresences of the two walkers, which we refer to as the inter-copresence time (ICT), by φ c (t). Then, the joint probability density of τ and n = 1 is given by

where S e (t) = e −λt and S d (t) = e −2Dt are the survival functions of φ e (t) and φ d (t), respectively, and S c (t) = ∞ t dt φ c (t ) and h c (t) = φ c (t)/S c (t) are the survival and hazard functions of φ c (t), respectively. By Laplace transforming p(τ, 1), we obtain

Similarly, the Laplace transform of the joint probability density of τ and n is given bŷ

Therefore, the IET distribution in the frequency domain is given bŷ

Using Eq. (6), we obtain the following expression of the coefficient of variation (CV), defined as the standard deviation divided by the mean, of τ (see Appendix A for the derivation):

whereφ c (0) andφ c (0) are the first and second derivatives ofφ c (s) evaluated at s = 0, respectively.

The asymptotic form of the distribution of first-passage time of random walks in networks with short relaxation time has an exponential tail [66] . With φ c (t) = αe −αt , which we refer to as the exponential ansatz, Eq. (6) is reduced top

The inverse Laplace transform of Eq. (8) is given by

where ∆ ≡ (λ + 2D + α) 2 − 4αλ > 0. Therefore, the distribution of IET is a mixture of two exponential distributions (see Appendix B for a proof). By substituting φ c (t) = αe −αt in Eq. (7), we obtain

FIG. 1. Schematic illustration of the dynamics of two walkers in a network with four subpopulations. Initially, at t = t0, the walkers are copresent in subpopulation 1. Then, events between them occur before t1. At t = t1, one walker moves from subpopulation 1 to subpopulation 4. At t = t2, the other walker moves from subpopulation 1 to subpopulation 2. At t = t3, the same walker moves to subpopulation 4 to be copresent with the other walker. Therefore, the first ICT is given by t3 − t1.

Although the walkers are copresent for t ∈ (t3, t4), no event occurs. At t = t4, one walker moves to subpopulation 3. At t = t5, the other walker also moves to subpopulation 3, producing a second ICT, which is equal to t5 − t4. Finally, an event occurs during their copresence in subpopulation 3, yielding an IET, denoted by τ * . In this example, this IET is followed by several much shorter IETs.

To mathematically derive the distribution of ICTs, φ c (t), we denote by s = (m, n) the state of the system in which one walker is in the mth subpopulation and the other walker is in the nth subpopulation. Because the walkers are indistinguishable, the system has N = N (N + 1)/2 states, i.e., (m, n), where 1 ≤ m ≤ n ≤ N .

We consider the first-passage time of a continuoustime Markov chain from a state (i, i), i.e., a copresence, to a state (j, j), i.e., the next copresence, where i, j ∈ {1, 2, . . . , N }. A similar problem of first-passage time to the next copresence of two random walkers appears in the theory of cooperation in evolutionary games on networks [67, 68] . We refer to the states in which the two walkers are in the same subpopulation as absorbing states and those in which the walkers are in different subpopulations as transient states. Let A = {(m, n); m = n} and B = {(m, n); m < n} be the set of absorbing and transient states, respectively. We denote the probability that the dynamics starting in state (i, i) at time 0 is in state s = (m, n), where m ≤ n, at time t by p i, s (t), or equivalently, p i,(m,n) (t). We set

At time 0, one of the walkers leaves subpopulation i such that the system leaves the absorbing state (i, i) and enters a transient state. Therefore, the distribution of the initial state of the Markov process is given by

where k out

A il is the out-degree of subpopulation i, and δ ij is the Kronecker delta.

The master equation for a transient state s = (m, n) is given by

where L ≡ 2DI −Q, matrix Q is the matrix of transitions rates between transient states, and I is the N T × N T identity matrix (see Appendix C for details). Equation (13) leads to

The master equation for the absorbing states is given by

where R is the matrix of transition rates from a transient state to an absorbing state (see Appendix C). By substituting Eq. (14) into Eq. (15), we obtain

The probability density function of the first-passage time from (i, i) to a state s, denoted by f i, s (t), satisfies

where f A i, s (t), or equivalently f A ij (t), is the probability density of the first-passage time from (i, i) to (j, j), because p s, s (t − t ) = 1 for any t ≥ t . Therefore, the probability density of the first-passage time to the absorbing states is given by

where

is the probability density of the first-passage time from (i, i) to (j, j).

We define the N × N T matrix P T (t) such that its ith row is given by

where P T 0 = P T (0). We are interested in the IET distribution in the equilibrium. In the equilibrium, the initial copresence state (i, i) should be the stationary probability of the discretetime random walk on A, where the transition probability T ij from (i, i) to (j, j) is given by

Using Eqs. (20) and (21), we obtain the N × N matrix

We obtain the stationary distribution q * for the absorbing states as

which we can numerically solve for arbitrary networks of subpopulations.

If Q is diagonalizable, we can write

where v j and w j are the right and left eigenvectors, respectively, associated with eigenvalue γ j of Q. By combining Eqs. (20) , (23) , and (24), we obtain the ICT distribution weighted by the stationary probability of the initial location of the two copresent walkers as

where c j = q * P T 0 v j w j 1 N T , α j = 2D − γ j > 0, and 1 N and 1 N T are the column vectors of size N and N T , respectively, in which all the elements are 1 [see Appendix D for a proof of Eq. (25) and that α j > 0]. Because

is a mixture of exponential distributions, including the case in which c j < 0.

By substituting the Laplace transform of Eq. (25) into Eq. (6) and Eq. (7), we obtain the IET distribution in the frequency domain and the CV aŝ

and

respectively. This solution depends on the network structure through φ c (t).

To compare CV ansatz with CV mixture we impose that the mean ICT is equal between the two. This condition combined with φ c (t) = αe −αt and Eq. (25) yields

Then, CV mixture is larger than CV ansatz because

2 , which is satisfied by the Sedrakyan's inequality (also known as Titu's lemma) [69] . Therefore, the exponential ansatz, which has led to the mixture of two exponential distributions for IETs, gives a lower bound in terms of the dispersion of IETs. Additionally, in the case of the network of two subpopulations connected to each other, Eq. (25) is reduced to the exponential ansatz (see Appendix D). Therefore, although the time-domain solution of the IET distribution is not available in general, the IETs are guaranteed to be distributed according to a distribution with a heavier tail than a mixture of two exponential distributions.

In fact, if the system starts in state (i, i) and reaches (j, j) to produce an ICT, it restarts from (j, j) to produce the next ICT. Therefore, Eq. (26) , which assumes the independence of different ICTs, is only approximate. Nevertheless, the following numerical simulations support that the difference between the approximate solution [i.e., Eqs. (26) and (27)] and the exact solution (see Appendix E) is negligible.

We simulated the model to validate our theory. We show the IET distribution produced by two walkers on the Barabási-Albert (BA) and Watt-Strogatz (WS) networks in Figs. 2(a) and 2(b), respectively. The figure shows that the IET distributions obtained from numerical simulations and ansatz have heavier tails than the exponential distribution whose mean IET is equal to that for the numerical simulations. The distribution obtained from the simulation decays more smoothly than the exponential ansatz.

To be more quantitative and general in terms of the parameter values and variety of networks, we compared the numerically and theoretically obtained CV of IET on

Survival function of IETs, S(τ ). We show the results for the direct numerical simulation of the model (simulation), the exponential ansatz [ansatz; Eq. (9)], and the single exponential distribution (exponential) whose mean is equal to that for the direct simulations. In Eq. (10), we set 10) and (27)], which only depends on D/λ. Therefore, we set λ = 1 without loss of generality and varied D.

The results shown in Fig. 3 indicate that the CV is large when D/λ is small for all four model networks and two empirical networks. Moreover, the CV obtained from our theory [Eq. (27) . Therefore, spatiality, or the large average path length between nodes, which is present in the last two networks but not in the other networks, may negatively affect the accuracy of the exponential ansatz.

Next, we investigated the effects of the network size on the IET. We set λ = 1 and D = 0.25, and computed the CV value for each of the six networks shown in Fig. 3 for different numbers of nodes, i.e., N = 50, 100, 200, 500, 1000, and 2000. We used the same parameter values as those used in Fig. 3 for each network model.

The CV values for the different network models and different N values are shown in Fig. 4 . We find that the results are roughly independent of N , except for the geographical threshold graph. For the geographical threshold graph, the CV is large at small N . Nevertheless, the CV value for the geographical threshold graph depends little on N for larger networks.

We have assumed simple random walks. Second-order random walks better approximate human and animal mo- [71] . (c) Geographical threshold (Geo. thld.) graph with the threshold value θ = 95, dimension equal to 2, Euclidean distance metric, and h(r) = r −2 [72] . (d) Lancichinetti-Fortunato-Radicchi (LFR) model with γ = 3, β = 1.5, µ = 0.2, average degree equal to 5, maximum degree equal to 50, and minimum community size equal to 10 [73] . (e) US airport network [74] . (f) Medellin intercity zone network [75] . In panels (a)-(d), the networks have N = 100 nodes. We set λ = 1.

bility [76, 77] . Our theory and the result that IETs obey a heavier-tailed distribution than the exponential distribution holds true for second-and higher-order random walks (see Appendix F).

More realistic metapopulation models incorporate the recurrent nature of human mobility where individuals iterate between a home subpopulation and a work subpopulation, which are different for different individuals [78] [79] [80] [81] [82] [83] . In this section, we analyze the IET distribution for such a mobility rule.

We emulate this situation by assigning to each walker a home subpopulation and allowing it to visit subpopulations adjacent to the home subpopulation. Then, each walker is confined to a star graph, of which the central subpopulation is the home subpopulation and the leaves are other locations such as work. Without loss of generality, we consider two random walkers on the respective star graphs to generate sequences of IETs.

We examine two cases. In the first case, the two walkers share the home subpopulation and all leaf subpopulations, as shown in Fig. 5(a) . In other words, the two walkers are confined to the same star graph, and events between them may occur in any subpopulation. Therefore, our theory directly applies.

In the second case, we assume that the two walkers have different home subpopulations and share some but not necessarily all the leaf subpopulations, as shown in Fig. 5(b) . This case mimics, for example, the situation in which the two individuals are coworkers living in different cities. Events may occur between the walkers only when they are copresent in any of the shared leaf subpopulations. To derive the states and the transition rates between the states, we distinguish the two walkers and the subpopulations. There are three types of subpopulations: the home subpopulations, which are unshared, the shared leaf subpopulations, which are work locations, for example, and the unshared leaf subpopulations, which represent other locations that one but not both walkers visits. In a network with N subpopulations, two subpopulations are home subpopulations, N s subpopulations are the shared leaf subpopulations, and walkers 1 and 2 have N 1 and N 2 unshared leaf subpopulations, respectively, such that N = 2 + N s + N 1 + N 2 .

We define the home subpopulation nodes by h 1 and h 2 , and the set of shared leaf subpopulations by W = {w 1 , w 2 , . . . , w Ns }. There are four types of transient states. In the first type of transient states, both walkers are in different shared leaf subpopulations. There are N s (N s − 1) such states. In the second type, walker 1 is in a shared leaf subpopulation and walker 2 is not There are N s (N 1 + 1) such states. In the third type, walker 2 is in a shared leaf subpopulation and walker 1 is not There are N s (N 2 + 1) such states. In the fourth type, both walkers are not in any of the shared leaf subpopulations. There are (N 1 + 1)(N 2 + 1) such states. Therefore, there are 

where χ W is the indicator function defined by

In an absorbing state, the two walkers are copresent in one of the shared leaf subpopulations. Therefore, there are N s absorbing states. Each element of the N T × N w matrix of transition rates from a transient state to an absorbing state, R s s , where s = (m 1 , m 2 ) ∈ B and s = (m 1 , m 2 ) ∈ A, is given by

The distribution of the initial state of the Markov process, is given by

Then, with Eqs. (30) , (32) , and (33) as inputs to our theory, the calculation steps to derive the distribution and CV of IET remain the same.

The CV values obtained from the simulation and the theory for the case in which the two walkers share the home and the case in which they have different homes are shown in Figs. 5(c) and 5(d), respectively. In both cases, the results are similar to those shown in Fig. 3 , i.e., the CV is large when D/λ is small, and the theory agrees with the simulation. The exponential ansatz solution yields slightly smaller values of CV than the numerical simulations, which we also observe in Fig. 3 .

We have shown that IET distributions are mixtures of exponential distributions for mobile individuals in metapopulation networks. The results hold true under mild conditions, i.e., for various structures of the metapopulation network, a work-home mobility rule, and higher-order random walks. Although a mixture of exponential distributions is technically not heavytailed, it often approximates heavy-tailed distributions reasonably well over scales [35, 37, 84, 85] . Therefore, the present results provide a compelling explanation of heavy-tailed IET distributions in human and animal contact data. Additional mechanisms such as circadian or weekly rhythms [33] and dynamics of individuals' internal states (e.g., high-activity versus low-activity states) [39] on top of mobility and metapopulation networks may make IET distributions more smooth and more powerlaw-like.

We assumed that two walkers in the same subpopulation meet each other according to a Poisson process. In other words, the waiting time until they have the next event obeys an exponential distribution if neither walker leaves the current subpopulation. In fact, this assumption is unnecessary for our results to hold true. One can assume other types of distributions for the mentioned waiting time to obtain qualitatively the same results. This is because the result that the IET distribution is a mixture of exponential distributions owes to the fact that the ICT distribution is a mixture of exponential distributions.

The present work bridges statistics of IETs and metapopulation models, which have been extensively but in most cases separately investigated topics in network science and related research fields. For example, heavy-tailed IET distributions suppress epidemic spreading in a majority of scenarios [7] [8] [9] [10] [11] [12] [13] [14] . The present results allow a new interpretation of these results. That is, different results of epidemic spreading associated with non-Poissonian IET statistics may owe in large part to whether the structure of the metapopulation network plays a significant role (which would yield heavier-tailed IET distributions) or the population is sufficiently wellmixed (which would yield an exponential IET distribution). The effects of IET statistics other than their heavy-tailed distributions [2, 3, 86] on epidemic and other dynamical processes may also be due to the underlying metapopulation network. Reexamining the role of IETs in contagion and other dynamical processes from the viewpoint of metapopulation networks warrants future work. The coefficient of variation (CV) is defined as

where · is the expectation. From Eq. (6), we obtain the first and second moments of τ as

and

(A3) In this section, we show that Eq. (9) is a mixture of two exponential distributions. We first rewrite Eq. (9) as

where

and

The two exponents in Eq. (B1) are positive because

the weights C 1 and C 2 are positive. Equations (B2) and (B3) also imply that

Therefore, Eq. (B1) is a mixture of two exponential distributions.

An element of the transition rate matrix W from a state s = (m, n) to a state s = (m , n ) is given by

The first term on the right-hand side of Eq. (C1) represents the transitions between transient states. We define Q as the N T × N T matrix of transition rates between transient states, i.e., W s s = Q s s , where s ∈ B and s ∈ B. Then,

The second term on the right-hand side of Eq. (C1) represents the transitions from a transient state to an absorbing state. We define R as the N T × N matrix of transition rates from a transient state to an absorbing state, i.e., W s s = R s s , where s ∈ B and s ∈ A. Then,

Note that

for any s ∈ B. Equation (C4) represents the fact that the system leaves any transient state at rate 2D owing to the movement of each walker, which occurs at rate D.

The master equation for a transient state s = (m, n), which is equivalent to Eq. (13) , is given by

In Eq. (C5), we used the fact that the sum of transition rates from any transient state to other states is equal to 2D, as shown in Eq. (C4).

To derive the ICT distribution, we first note that Eq. (C4) leads to

Then, by combining Eqs. (20) , (23) , and (D1), we obtain

By substituting Eq. (24) into Eq. (D2), we find

In our model, all transient states are reachable from any transient state. Therefore, Q is irreducible and, by definition, non-negative. Hence, by the Perron-Frobenius theorem, the spectral radius of Q, denoted by ρ(Q), satisfies ρ(Q) ≤ max s s ∈B Q s s . Using Eq. (C4), we obtain ρ(Q) ≤ 2D, which implies that 2D − γ j ≥ 0, for all j.

In the particular case of two subpopulations connected to each other (N = 2), the number of transient states is N T = 1. Therefore, Eq. (D3) becomes

which is equivalent to the exponential ansatz. In this case, there is one transient state and two absorbing states. Therefore, Q = 0, R = D D , L −1 = 1/2D, and 1 N T = 1. From Eq. (22), we obtain

such that q * = 1/2 1/2 . Then, we obtain α 1 as

Therefore, for the case of two subpopulations, the ICT distribution is given by

Appendix E: Exact solution for the interevent time distribution

In this section, we derive the exact solution for the IET distribution. The probability density with which the two walkers are copresent in subpopulation j after time t given that the last copresence terminated at time 0 in subpopulation i is given by f A ij (t), which is given via Eq. (20) . In other words, f A ij (t) is the probability density of an ICT, t, from (i, i) to (j, j). Note that f A ij (t) = φ c (j, t|i) and the normalization is given by N j=1

∞ 0 dt f A ij (t) = 1. We denote by p ij (τ, n) the joint distribution of an IET, τ , and the number of copresences, n, such that the event has occurred in subpopulation i and the next event occurs in subpopulation j. The normalization is given by N j=0

∞ n=0 ∞ 0 dτ p ij (τ, n) = 1. We derive the Laplace transform of p ij (τ, n) by extending Eq. (5) as follows:

We define the N × N matrixΦ(s) by [Φ(s)] ij = φ c (j, t|i). Using Eq. (20), we obtain

We also define the N × N matrixΠ(s, n) by [Π(s, n)] ij =p ij (s, n). Using Eq. (E1), we obtain

(E3)

We denote byĝ(s) the stationary distribution in the frequency domain of IET weighted by the stationary probability of the initial location of the two copresent walkers. By combining Eqs. (E2) and (E3), we obtain

where q * is given by Eq. (23). With Eq. (E4), the first and second moments of IET are given by

and

respectively, wherê

By substituting Eqs. (E5) and (E6) into Eq. (A1), we obtain

Appendix F: Higher-order random walks

In this section, we show that our theory holds true when the individuals move according to higher-order random walks. We focus on second-order random walks and then explain how our theory generalizes to higher orders.

In a second-order random walk, the probability with which a walker visits the next node depends on its current and last visited nodes [76, 77] . Therefore, the state of a single walker is defined by a pair of nodes (m − , m), where m is the currently visited node, and m − is the node that the walker visited just before arriving in m. It should be noted that m = m − and that we distinguish between (m − , m) and (m, m − ). In other words, the state of each walker is specified by a directed edge. For example, the sequence of a walker's positions from node i to j and then to k is given by (i, j) → (j, k). Therefore, in a second-order random walk, we can regard the movement of the walkers as a first-order random walk from directed edge (i, j) to directed edge (j, k) instead of between nodes. We assume in the following text that there are N subpopulations and M undirected edges. Then, a walker moves among 2M directed edges.

The state of the system of two second-order random walkers is defined by a pair of directed edges, one for each walker. For example, if one walker is currently in subpopulation m and the other walker in subpopulation n, we denote the state of the system by ((m − , m), (n − , n)). Because the walkers are indistinguishable the system has N = 2M (2M + 1)/2 = M (2M + 1) states in total.

We denote the next state of the system by ((m − , m ), (n − , n )) and suppose that the walker at subpopulation m has moved to subpopulation m new . Therefore, the new state ((m − , m ), (n − , n )) is equal to ((m, m new ), (n − , n)). If m new = n, then the next state (i.e., ((m − , m ), (n − , n ))) is an absorbing state. Otherwise, it is a transient state.

The number of absorbing states with which the two walkers meet at node i is given by k i + ki 2 = k i (k i +1)/2, where k i is the degree of node i. Note that ((j, i), (j, i)), where j is a neighbor of i, is also a valid absorbing state. Then, the system has

i absorbing states. The number of transient states is given by

i . Because we have identified all the transient states, absorbing states, and state transition rules, one is able to define the N T × N T matrix Q of transition rates between transient states and the N T × N A matrix R of transition rates from a transient state to an absorbing state. Then, the theory that follows is the same as that developed in Section III C and in Appendix E.

It is straightforward to extend the same procedure to the case of directed metapopulation networks and higher-order (i.e., third-order or higher) random walks. For a third-order random walk, for instance, the state of the system is described by a pair of triples, i.e., ((m −2 , m −1 , m), (n −2 , n −1 , n)), where m is the subpopulation that the first walker currently visits, m −1 is the node that the same walker visited just before m, and m −2 the subpopulation that the walker visited just before m −1 . The definitions are analogous for (n −2 , n −1 , n).

As an example of a second-order random walk, we simulated the non-backtracking random walk. By definition, a non-backtracking random walker on an undirected unweighted network that has moved from m − to m, moves to any of the neighbors except m − with the equal probability in the next move [87, 88] . The CV of IET for two non-backtracking random walkers on various networks is shown in Fig. 6 . The results are similar to those for the simple random walk shown in Fig. 3 . In other words, the CV is substantially larger than 1 and large when D/λ is small. In general, the CV values for the non-backtracking random walk are somewhat smaller than those for the simple random walk. 

Temporal networks

A Guide to Temporal Networks

Bursty Human Dynamics

Temporal Network Theory

The origin of bursts and heavy tails in human dynamics

Modeling bursts and heavy tails in human dynamics

Spreading dynamics following bursty human activity patterns

Small but slow world: How network topology and burstiness slow down spreading

Simulated epidemics in an empirical spatiotemporal network of 50,185 sexual contacts

Dynamical strength of social ties in information spreading

Predicting and controlling infectious disease epidemics using temporal networks, F1000Prime Rep

Analytically solvable model of spreading dynamics with non-Poissonian processes

Epidemic processes in complex networks

Temporal Network Epidemiology

Human comment dynamics in on-line social systems

Voter model with non-Poissonian interevent intervals

Update rules and interevent time distributions: Slow ordering versus no ordering in the voter model

Dynamics of social balance under temporal interaction

Evolution of cooperation on temporal networks

Threshold model of cascades in empirical temporal networks

Bursty communication patterns facilitate spreading in a thresholdbased epidemic dynamics

Effects of temporal correlations on cascades: Threshold models on temporal networks

Dynamics of cascades on burstiness-controlled temporal networks

Generalized master equations for non-Poisson dynamics on networks

Random walks on temporal networks

Steady state and mean recurrence time for random walks on stochastic temporal networks

Random walks and diffusion on networks

Exact results for the Barabási model of human dynamics

Biased diffusion and universality in model queues

Priority queues with bursty arrivals of incoming tasks

Impact of interactions on human dynamics

Time-varying priority queuing models for human dynamics

A Poissonian explanation for heavy tails in e-mail communication

On universality in human correspondence activity

Selfexciting point process modeling of conversation event sequences

Small inter-event times govern epidemic spreading on networks

Long-tailed distributions of inter-event times as mixtures of exponential distributions

Two-state markov-chain Poisson nature of individual cellphone call statistics

Generative models of simultaneously heavy-tailed distributions of interevent times on nodes and edges

Predictability of conversation partners

On the dynamics of human proximity for data diffusion in ad-hoc networks

Fingerprinting temporal networks of close-range human proximity

Temporal networks of face-toface human interactions

Activity clocks: Spreading dynamics on temporal networks of human contact

Contact patterns among high school students

How memory generates heterogeneous dynamics in temporal networks

Can co-location be used as a proxy for face-to-face contacts?

Modeling human dynamics of face-to-face interaction networks

Model reproduces individual, group and collective dynamics of human contact networks

Modelling temporal networks of human face-to-face contacts with public activity and individual reachability

Emergence of metapopulations and echo chambers in mobile agents

Similarity forces and recurrent components in human face-to-face interaction networks

Reaction-diffusion processes and metapopulation models in heterogeneous networks

An immunization model for a heterogeneous population

Spatial heterogeneity and the design of immunization programs

Spatial heterogeneity in epidemic models

(Meta) population dynamics of infectious diseases

Cities and villages: Infection hierarchies in a measles metapopulation

Metapopulation dynamics

Metapopulation Biology

Genetics and Evolution of Metapopulations

Invasion threshold in heterogeneous metapopulation networks

Epidemic modeling in metapopulation systems with heterogeneous coupling pattern: Theory and simulations

Dynamical Processes on Complex Networks

The effect of travel restrictions on the spread of the 2019 novel coronavirus (COVID-19) outbreak

Asymptotic analysis of first passage time in complex networks

Evolutionary dynamics on any population structure

Evolution of cooperation on large networks with community structure

Algebraic Inequalities

Emergence of scaling in random networks

Collective dynamics of 'small-world' networks

Geographical threshold graphs with small-world and scale-free properties

Benchmark graphs for testing community detection algorithms

Rich do not rise early: Spatio-temporal patterns in the mobility networks of different socio-economic classes

Causality-driven slow-down and speed-up of diffusion in non-markovian temporal networks

Memory in network flows and its effects on spreading dynamics and community detection

Phase transitions in contagion processes mediated by recurrent mobility patterns

Natural human mobility patterns and spatial spread of infectious diseases

Human mobility and time spent at destination: Impact on spatial epidemic spreading

Critical regimes driven by recurrent mobility patterns of reaction-diffusion processes in networks

Epidemic spreading in localized environments with recurrent mobility patterns

Impact of temporal scales and recurrent mobility patterns on the unfolding of epidemics

Fitting mixtures of exponentials to long-tail distributions to analyze network performance models

Latent geometry and dynamics of proximity networks

Randomized reference models for temporal networks

Nonbacktracking random walks mix faster

Non-backtracking random walk