Analytically Solvable Model of Spreading Dynamics with Non-Poissonian Processes

Non-Poissonian bursty processes are ubiquitous in natural and social phenomena, yet little is known about their effects on the large-scale spreading dynamics. In order to characterize these effects we devise an analytically solvable model of Susceptible-Infected (SI) spreading dynamics in infinite systems for arbitrary inter-event time distributions and for the whole time range. Our model is stationary from the beginning, and the role of lower bound of inter-event times is explicitly considered. The exact solution shows that for early and intermediate times the burstiness accelerates the spreading as compared to a Poisson-like process with the same mean and same lower bound of inter-event times. Such behavior is opposite for late time dynamics in finite systems, where the power-law distribution of inter-event times results in a slower and algebraic convergence to fully infected state in contrast to the exponential decay of the Poisson-like process. We also provide an intuitive argument for the exponent characterizing algebraic convergence.


I. INTRODUCTION
Events of the dynamical processes of various complex systems are often not distributed homogeneously in time but have intermittent or bursty character. This is ubiquitously witnessed in processes of nature like earthquake statistics [1], solar flare [2], and firing of neurons [3], but also in social processes like financial interactions [4] and human communication activities [5]. In these examples, the distribution of inter-event times follows a power-law behavior [1][2][3][4][5], in contrast to the homogeneous Poissonian processes showing exponential distribution.
Dynamical processes of complex systems can be considered to take place on a network formed by pairwise interactions between the constituents of the system [6,7]. In the recently developed approach of temporal networks [8], a link between two nodes is considered to exist only at the moment of interaction. One of the most interesting dynamical processes on networks is the spreading [9][10][11][12][13] that takes place on temporal networks, and the statistics of events strongly influences its most important feature, namely, the speed of propagation. This feature is of pivotal interest and importance, e.g., for halting epidemic outbreaks or promoting diffusion of innovations.
Recently much effort has been devoted to clarifying how burstiness of events influences the spreading speed, partly by using empirical data analysis [14][15][16][17][18][19] and partly by model calculations [8,14,[20][21][22]. In those studies, the bursty character of an event sequence was found to slow down the late-time dynamics of spreading, evidenced also by a heavy tail in the inter-event time distribution. However, for the early-time dynamics, conflicting results have been reported [23]. In studies by Vazquez et al. [14] and Karsai et al. [15], the burstiness is found to slow down spreading, while other works point towards the opposite direction [16,18,21]. This calls attention to the importance of small inter-event times or to the role of lower bounds of inter-event times. The effect of lower bounds has been largely ignored, although it is present in empirical phenomena and it is important for understanding the early-time behavior of reference systems.
Model studies usually aim to reproduce some empirical observations to uncover the main mechanisms of the real underlying processes. Here, we take the perspective of constructing a tractable and analytically solvable model, where the effects of different parts of the inter-event time distribution can be studied explicitly and understood in detail. The model we consider is without correlations found in realistic data sets, except for the correlation due to the inter-event time distribution. In this way, we hope to provide a reference system, which can serve as a starting point for later studies.

II. MODEL
In order to model bursty spreading phenomena, we study deterministic susceptible-infected (SI) dynamics taking place on a temporal network of infinite size. Each node remains inactive for an inter-event time, denoted by l, before becoming instantaneously active, and then it turns inactive for another inter-event time period, and so on. The inter-event time distribution P 0 ðlÞ is assumed to be the same for all nodes, implying a homogeneous population. The activation pattern of a node is independent of whether it is susceptible (S) or infected (I). Whenever any infected node becomes active, it chooses a random node and infects it if the chosen node is susceptible (see Fig. 1). Here, the probability of choosing a susceptible node is 1 in the infinite-size system, as the dynamics starts from a single infected node. The newly infected inactive node should wait a residual time, denoted by w, before it becomes active again. The distribution of w is derived from P 0 ðlÞ as P 1 ðwÞ ¼ 1 μ R ∞ w P 0 ðlÞdl, with μ denoting the mean of l. Thus, the dynamics is stationary from the beginning, as it is independent of the nodes being initially active or inactive. Otherwise, the early stage of spreading dynamics could be sensitive to the variation of the initial distribution of active or inactive nodes.
In our model, the dynamics can be interpreted as occurring on a temporal network, in the sense that any pairwise interaction between nodes, defining a link, is instantaneous and annealed. Such links can be interpreted as directed, as the inter-event time distribution is considered only for outgoing events of infecting nodes. The spreading dynamics on this temporal network can be related to a class of Bellman-Harris branching processes [16,20,24]. It should be noted that temporal inhomogeneities have been considered in a model study by Perra et al. [12], although they took a different approach than ours by using inhomogeneous activities of nodes.
We investigate the spreading dynamics, starting from one infected and active node at time t ¼ 0. The number of infected nodes at a later time t, denoted by I 0 ðtÞ, remains 1 for the inter-event time l given to the initially infected node. After the first infection at t ¼ l, I 0 ðtÞ can be written as the sum of two numbers: One is for the infecting node and its subsequent infected nodes, which can be denoted by an independent and identical copy of I 0 but starting at t ¼ l, thus as I 0 0 ðt − lÞ. The other is for the newly infected node and its subsequent infected nodes, similarly denoted by Here, I 0 1 is an independent and identical copy of I 1 , i.e., the number of infected nodes starting from one infected and inactive node. Thus, we get (1) as depicted in Fig. 1(a). Since the newly infected node must wait a residual time w, as in Fig. 1(b), the number of infected nodes starting from one infected and inactive node can be written as where the I 00 s are independent and identical copies of I. The generating function for I 0 ðtÞ is defined as F 0 ðz; tÞ ¼ P k≥0 Pr½I 0 ðtÞ ¼ kz k , and we get Here, F 1 ðz; tÞ is the generating function defined for I 1 ðtÞ. By taking the expectation over l with P 0 ðlÞ, one obtains Then, the average number of I 0 ðtÞ is calculated as where n 1 ðtÞ ≡ hI 1 ðtÞi. Taking the Laplace transform, we getñ which straightforwardly leads tõ where we have used the relationP 1 ðsÞ ¼ 1 μs ½1 −P 0 ðsÞ. Then, n 0 ðtÞ can be calculated by taking the inverse Laplace transform ofñ 0 ðsÞ analytically or numerically if necessary. Note that this solution has been obtained for arbitrary interevent time distributions and for the whole time range, enabling us to evaluate the effect of burstiness on spreading at any stage of dynamics. In contrast to this case of infinite system size, the late-time behavior of finite systems cannot be investigated analytically.
As for the non-Poissonian bursty processes, they are often characterized by broad inter-event time distributions, such as Gamma and log-normal distributions [16] and power-law distributions with an exponential cutoff [14,21]. Since these distributions have a zero lower bound for interevent times, the effect of the lower bound on the early stage of spreading dynamics has been ignored, despite the importance of the finite lower bound in empirical phenomena. In order to systematically investigate the effect of the lower bound as well as the heavy tails of inter-event times, we consider the shifted power-law distribution with an exponential cutoff: where Γ is the upper incomplete Gamma function and θ is the Heaviside step function. l 0 and l c denote the lower bound and the exponential cutoff of l, respectively. In the case with l c → ∞, the value of the power-law exponent α should be larger than 2 to guarantee finite μ, i.e., the mean of l. The mean is related to other parameters as follows: where x ≡ l 0 μ and y ≡ l 0 l c . Here, 0 ≤ x ≤ 1 because the mean cannot be smaller than the lower bound. When y ¼ 0, the relation reduces to x ¼ α−2 α−1 . Note that setting the power-law exponent α ¼ 0 reduces the distribution to the shifted Poissonian case.

A. Poissonian processes
As the simplest case, the Poissonian process with P 0 ðlÞ ¼ μ −1 e −l=μ results in the solution n 0 ðtÞ ¼ e t=μ : Generally, we consider the shifted Poissonian process by setting α ¼ 0 in Eq. (10), leading to Here, we have used l c ¼ μ − l 0 from Eq. (11). Then, we getñ For the early-time dynamics, by assuming that s ≫ 1, we obtainñ which results in the following solution: As for Eq. (17), the exponential functions for t ≪ μ have been expanded. The lower bound l 0 delays the first branching, while at times later than t ¼ l 0 , the spreading speeds up. Let us define a dimensionless spreading rate at the moment of the first branching as follows: Next, we study the late-time dynamics, where the late time for infinite-size systems corresponds to the intermediate time for finite-size systems. Since it is evident that n 0 ðtÞ ∼ e t=μ for large t, we characterize the asymptotic behavior by defining the coefficient of the leading exponential term as Here, we have used the final value theorem [25] lim t→∞ fðtÞ ¼ lim s→0 sfðsÞ; with fðtÞ ¼ n 0 ðtÞe −t=μ . By plugging Eq. (14) into Eq. (21), we obtain This result implies that the finite lower bound suppresses the late-time spreading dynamics.

B. Non-Poissonian processes
Now we consider the general form of P 0 ðlÞ in Eq. (10). The Laplace transform of P 0 ðlÞ is obtained as follows: To investigate the early-time dynamics of n 0 ðtÞ, we consider the case with s ≫ 1. By expanding the incomplete gamma function, we obtaiñ Γð1−α;yÞ , leading to The spreading rate at t ¼ l þ 0 , i.e., C 0 , is obtained as Note that x, y, and α are not independent by means of Eq. (11) and that C 0 ðx; y; 0Þ ¼ 1 1−x ¼ C P 0 ðxÞ. For the comparison to shifted Poissonian processes, we define the ratio of spreading rates as r 0 ðx; y; αÞ ≡ C 0 ðx; y; αÞ C P 0 ðxÞ ; which turns out to be exactly the same as the ratio of This indicates that the probability of having l ¼ l þ 0 determines the early-time spreading dynamics. In addition to r 0 ðx; y; 0Þ ¼ 1, by definition, it is found that r 0 > 1 for α > 0 and r 0 < 1 for α < 0 [see Fig. 2(a)]. Provided that α > 0, one can conclude that the non-Poissonian bursty activity always accelerates the early-time spreading dynamics as compared to the shifted Poissonian case with the same mean μ and the same lower bound l 0 .
For the late-time dynamics, we focus on the asymptotic behavior characterized by C ∞ in Eq. (20). Similarly to the Poissonian case, one gets the following general result: Similarly to r 0 , we define the ratio r ∞ ðx; y; αÞ ≡ C ∞ ðx; y; αÞ C P ∞ ðxÞ : In addition to r ∞ ðx; y; 0Þ ¼ 1, by definition, it is found that r ∞ > 1 for α > 0 and r ∞ < 1 for α < 0 [see Fig. 2(b)]. This implies that the non-Poissonian bursty activity accelerates the late-time dynamics as compared to the corresponding Poissonian processes. Our analysis is confirmed by the numerical simulations as depicted in Fig. 3(a). The correction term to the exponential growth for latetime dynamics is obtained for the case of y ¼ 0 such that Both B and λ can be analytically obtained by taking the Laplace transform of the above equation for s > μ −1 . By defining ϵ ≡ s − μ −1 > 0, Eq. (33) is transformed to By expanding both sides up to the order of ϵ 2 for small ϵ and comparing the coefficients, one can get λ and B as functions of x and α. In the limit of x → 0 (α → 2), we obtain λ → μ −1 and B → 1, resulting in  29) and (b) the ratio of asymptotic coefficients r ∞ ðx; y; αÞ in Eq. (32). We note that x ¼ l 0 μ and y ¼ l 0 l c , with the lower bound l 0 , the mean μ, and the cutoff l c of the inter-event time distribution in Eq. (10). We plot r 0 ðx; yÞ and r ∞ ðx; yÞ for various values of y and, in the insets, r 0 ðy; αÞ and r ∞ ðy; αÞ for various values of α. By definition, the ratios have the value of 1 for α ¼ 0. In the limits of x → 1 and/or y → ∞, one gets P 0 ðlÞ → δðl − μÞ; thus, r 0 and r ∞ have the value of 1.
The constant term implies the existence of residual times that are effectively infinite because of the non-normalizability of P 1 ðwÞ, whose tail is characterized by the exponent α − 1.

C. Finite-size effects
Finally, we consider the effect of the finite system size N on the spreading dynamics. Whenever an infected node becomes active at time t, it chooses a random node and infects it if the chosen node is susceptible. In other words, the infection occurs with probability N−n 0 ðtÞ N−1 ≈ 1 − n 0 ðtÞ N , but otherwise it does not occur. Since the exact solution could not be obtained, we perform numerical simulations to obtain the spreading dynamics shown in Figs. 3(b), 3(c), 3(d), and 3(e). The early-and intermediate-time dynamics are consistent with the early-and late-time dynamics of the infinite system, respectively. For the late-time dynamics, the non-Poissonian bursty activity results in a slower algebraic convergence to the fully infected state, i.e., ∼t −β with β ¼ α − 1 for the power-law case, in contrast to the exponential decay of the Poisson-like process, i.e., ∼e −t .
We provide an intuitive argument for the relation β ¼ α − 1 in the case of power-law inter-event time distributions. While the average fraction of susceptible nodes decays algebraically, the fraction of susceptible nodes for each run remains close to 1 and then suddenly decays exponentially, as shown in Fig. 3(e). The period staying close to 1 must be governed mostly by the inter-event time initially given to the first infected node. Therefore, the average fraction of susceptible nodes can be obtained as the fraction of runs that did not reach the fully infected state at time t. Such a fraction of runs is equal to the probability of having l > t, which is proportional to t −ðα−1Þ , leading to β ¼ α − 1. Considering the dominant role of inter-event time given to the first infected node, this argument should be valid for the SI dynamics with a power-law inter-event time distribution on a broad class of networks.

III. CONCLUSIONS
We have introduced an analytically solvable model for studying the effect of non-Poissonian bursty inter-event time distributions on the SI spreading dynamics. Our model is devised to be stationary from the beginning. For this, we make a realistic assumption that for each infection event the infecting node should wait another inter-event time and the newly infected node should wait a residual time before they become active again. With this assumption, we can obtain the analytic solution of spreading dynamics in infinite systems for arbitrary inter-event time distributions, but more importantly for the whole time range. With our analytic solution, the role of the lower bound of interevent times has been exactly compared for Poisson-like and non-Poissonian processes. We also note that, as done in the case of null models, randomizing or shuffling the event timings to destroy temporal correlations can eliminate the lower bound of inter-event times. Hence, for systematic comparison between the original situation and the null model, one needs to employ a shuffling method that conserves the lower bound.
Let us next discuss apparently conflicting results for the early-stage spreading dynamics presented in Refs. [15,18]. The early-stage dynamics is mainly driven by small interevent times, which generally leads to a faster spreading for non-Poissonian cases than for Poisson-like cases. This is the case when any infected node finds a susceptible node without a topological limit. However, it is well known that the mobile-phone-call network (MCN) in Ref. [15] has the community structure accompanying the bottleneck effect due to weak links between communities [26]. The large inter-event times associated with such weak links are the In (a,b), we set μ ¼ 1 and l 0 ≈ 0.106 for all cases, from which we get α ≈ 2.118 for l c → ∞ (power law) and α ¼ 1.5 for l c ¼ 100l 0 (power law with cutoff). Each curve is averaged over up to 1.5 × 10 6 runs. In (c), we plot the fractions of susceptible nodes, 1 − n 0 ðtÞ N , using the same results as in (b), as well as for various values of α in the power-law case with μ ¼ 1 (thin reddish curves). The Poisson-like case is fitted with e −t , while the power-law cases are fitted with t −β . The estimated β as a function of α is presented in (d), comparable with β ¼ α − 1. For understanding this scaling relation, in (e) we plot the individual runs with different initial inter-event times (from blue to red curves) with their average (black solid curve) in the case with α ≈ 2.118. main reason for slowing down the spreading on the MCN, while the role of small inter-event times coupled with local topological structure is still important for spreading within communities. This leads us to study more realistic models as a future work. On the other hand, the spreading is enhanced by the burstiness on the sexual network [18]. This might be because the sexual network has different community and/or temporal structures from the MCN, so the spreading on the sexual network can be better understood by our model to some extent.
As a followup, our model can be extended to incorporate a number of other complex situations, such as susceptibleinfected-recovered spreading dynamics and cascading phenomena. Our results should be of interest beyond the community of network scientists because non-Poissonian processes are ubiquitous and yet little is known about their impact on the large-scale dynamics.