The nature of epidemic criticality in temporal networks

Analytical studies of network epidemiology almost exclusively focus on the extreme situations where the time scales of network dynamics are well separated (longer or shorter) from that of epidemic propagation. In realistic scenarios, however, these time scales could be similar, which has profound implications for epidemic modeling (e.g., one can no longer reduce the dimensionality of epidemic models). We build a theory for the critical behavior of susceptible-infected-susceptible (SIS) epidemics in the vicinity of the critical threshold on the activity-driven model of temporal networks. We find that the persistence of links in the network leads to increasing recovery rates reducing the threshold. Dynamic correlations (coming from being close to infected nodes increases the likelihood of infection) drive the threshold in the opposite direction. These two counteracting effects make epidemic criticality in temporal networks a remarkably complex phenomenon.

The quest for epidemic models that are both realistic, analytically tractable, and realistic is a fundamental challenge to theoretical research.This is never as difficult as when one cannot rely on ignoring a structure, by assuming it constant.When building an analytical theory, how fast the network changes and how fast the epidemic spreads are two dynamics that one would like to reduce to one.From a medical point of view, that might not work-an incipient outbreak might sweep over the population in a matter of weeks, the same time scale of updates to the friendship network [1,4].This is the motivation for temporal network epidemiology, where the networks do not necessarily change slower or faster than the disease propagation [13,14].
There are three main philosophies of how to go beyond static network epidemiology to include the temporality of edges.The first approach temporal contacts as possible contagion events [8].This approach often relies on empirical data of the contacts and analyzes them by computer simulations.The second approach, adaptive networks, models the contacts without real-world data but assumptions about how they change given the state of the disease [6].In this paper, we will take a third approach which is also purely model-based, trying to emulate dynamic contact structures that real epidemics spread on, while being tractable for analytical calculations.Typically these models-such as blinking networks [7] and activity-driven networks [12,18]-assume an underlying static network over which they generate active contacts.
The original activity-driven network model was proposed by Perra et al. [18], as follows: (i) The network is initialized to be completely disconnected and an activity a i = ηx i is assigned to each node i based on the given activity distribution F(x); (ii) At each discrete time t, each node i becomes active with probability a i ∆t and randomly creates m undirected links with other nodes.(iii) The epidemic dynamics take place over the instantaneous network; (iv) At the next time step t + ∆t, all the links are removed, and the process resumes from step (ii).
In static network epidemiology, the canonical models (SIR and SIS) are effectively governed by one parameter β/µ-the ratio of infection to recovery rate.Therefore, most studies of activity-driven networks also focused on this ratio [5, 16-19, 21, 24, 25, 28].This approach, effectively studying the model in the limit of fast network dynamics, goes against the idea of temporal networks as the modeling framework for intermediate time scales.In this limit, the β/µ-based mean-field approach accurately describes the criticality of compartmental models on activity-driven networks.However, reality can have slower network dynamics, and then the mean-field analysis will fail.
More precisely: The thresholds of both annealed (the limit of rapidly changing) and static (the limit of slowly changing) networks depend only on β/µ, but their actual critical values (β/µ) c are different because of dynamical correlations in the static networks [3].(I.e., infected nodes tend to group together.)Since the activity-driven network can interpolate between these extremes, by keeping one of β or µ fixed, we can understand that tuning only the other must change the threshold (contradicting the assumption that β/µ fully describes the criticality of epidemics on the activity-driven model).
In this paper, we employ a continuous time description [11,29] of the coevolution of the network and SIS process.We use simulations to see the effects of both network correlations (reflecting the persistence of social contacts) and dynamic correlations caused by the inertia of the disease propagation.We complement these simulations with theoretical calculations that, for the first time in the literature, include arXiv:2311.12175v1[physics.soc-ph]20 Nov 2023 the above-mentioned network correlations but ignore dynamic correlations.The comparison of these two types of results allows us to see the effects of the two types of correlations.
Next, we turn to a technical definition of the model-a continuous time description of the coevolution of the network and SIS dynamics.For the evolution of the network, each inactive (U) individual i is activated with a rate a i = ηx i , where x i is drawn from a probability distribution F(x).When an individual gets activated, it randomly selects m individuals to generate links to.Active (A) individuals become unactivated with a rate b, and then remove all links to themselves.Here, the use of U instead of I for inactive individuals is to avoid confusion with the concept of the epidemic model.For the evolution of disease, we employ the classical SIS model.Susceptible (S) individuals become infected via contacts with infected (I) individuals at rate β times the number of S-I links.Infected individuals recover to susceptible with rate µ.The detailed simulation procedure see Appendix A.
We start by analyzing the network evolution, which is independent of the SIS model.Define ρ A a as the fraction of active individuals with activity rate a, and ρ A = ρ A a da as the total fraction of active individuals.The dynamic equation of the fraction of active individuals of class a in the network can be written as Here, the first term of Eq. ( 1) represents the spontaneous creation of active individuals, while the second term represents the spontaneous annihilation.When the evolution of the network reaches its steady state, i.e. dρ Aa dt = 0, we have where a a+b = ηxF(x) ηx+b dx.We define ⟨k A,a ⟩ as the average degrees of active individuals in class a. Considering that the links created per unit time in the network system should be equal to the disconnected ones, we have Combining Eqs. ( 2) and ( 3), we obtain where ⟨k A ⟩ is the average degrees of active individuals.Equation (4) also implies that the average degrees of inactive individuals in different classes are the same, and the average degree of inactive individuals is represented by ⟨k U ⟩.
Similarly, the number of links between active and inactive individuals is stable when the system reaches its steady state.Since aN(1 − ρ A a )da individuals are activated per unit time, the number of new links to the inactive individuals is m(1 − ρ A ) aN(1 − ρ A a )da.On the one hand, links are disconnected when the active node of the link deactivates.On the other hand, when an inactive node in the link is activated, the class of the link changes.Therefore, the amount of link reduction between active and inactive individuals per unit time is ⟨k U ⟩ N(1 − ρ A a )(a + b)da.Note that the number of links between active and inactive individuals is equal to ⟨k U ⟩ N(1 − ρ A a )da, because there are no links between inactive and inactive individuals.According to the relation Next, we consider the coevolution of networks and epidemics.Since individuals are infected through links, the evolution of the disease cannot be separated from the evolution of the network.The individuals are classified as XY, where X ∈ {U, A} and Y ∈ {S , I}.We introduce the notation ρ XY a for the fraction of individuals with state XY and class a; is the probability that an individual with state XY and class a is connected to an individual with state X ′ Y ′ and class a ′ .Note that ϕ US a ,US a ′ , ϕ UI a ,UI a ′ , and ϕ US a ,UI a ′ would be impossible.
From the above, we can state the master equations as follows, Since only active individuals can make connections and inactive individuals cannot, we can approximate ϕ XY a ,X ′ Y ′ a ′ in Eq. (6a) in the form of ρ XY a using the mean field anzats, i.e.
Note that the approximation of Eq. (7a) means that we ignore the correlation of dynamic.By inserting Eq. (7a) into Eq.(6a), and then integrating both sides of Eq. (6a) with respect to a, we get where ρ UI = ρ UI a da and ρ AI = ρ AI a da.See Appendix C for detailed derivations.Now we perform a linear stability analysis for Eq. ( 8) around the fixed point ρ AI = ρ UI = 0, and after ignoring the higher order terms, the Jacobian matrix of Eq. ( 8) can be written as When the largest eigenvalue of the Jacobian matrix is null, we can obtain the epidemic threshold When the probability distribution F(x) satisfies the Delta distribution, i.e., F(x) = δ(x − x 0 ), G and H in Eq. ( 10) can be simplified to G = (ηx 0 +b)(µ+ηx 0 +b) bmµ and H = bµ ηx 0 (µ+ηx 0 +b) .It is clear that the epidemic threshold of Eq. ( 10) is dependent on µ on the activity-driven network, which is an important difference from annealed and static limits.Specifically, as µ ≫ a and µ ≫ b, the activity-driven network degenerates into the static network, and the threshold of Eq. ( 10) can be rewritten as Besides, as µ ≪ a and µ ≪ b, the activity-driven network will degenerate into the annealed network, for which the threshold of Eq. ( 10) approximate by using Taylor expansion.
Near the epidemic threshold, the order parameter fluctuates greatly, leading to a critical slowing.The latter means that the relaxation time of the order parameter tends to infinity in the thermodynamic limit.In general, the internal causes of divergence in relaxation time and divergence in thermodynamic susceptibility are the same.To estimate the threshold of the SIS model on simulation, a lifetime-based method is proposed in Ref. [2] on static networks, and Ref. [19] use this method on activity driven networks.The lifetime L is defined as the time that passes before the disease either goes extinct or spreads to a finite fraction C of the network, where C is the fraction of distinct nodes ever infected during the simulation.When β/µ is much less than its critical value (β/µ) c , the disease will quickly die out.If β/µ ≫ (β/µ) c , the disease rapidly reaches a steady state.For β/µ ∼ (β/µ) c , we can observe a longer lifetime due to critical slowing down.In this case, the opposing mechanisms of infection and spontaneous recovery almost balance out, making spreading dynamics slow.In the thermodynamic limit, the average lifetime diverges at (β/µ) c both from below and from above.For finite systems, we use the maximum value of L to estimate the epidemic threshold.
The simulation results of Fig. 1, show an interesting nonmonotonic increase of the threshold with µ.Specifically, the Eq. ( 10 5).Scatters are simulation results with an accuracy of 0.01.Lines are the theoretical estimations obtained from Eq. ( 10) that consider only the network correlations.
curves show both local maxima and minima, whereas the analytical results are strictly decreasing.Meanwhile, the threshold tends to have different saturation values in the static or the annealed limits.The simulation results will increase with the increase of µ when the analytical results tend to be saturated.Based on the above analysis, one can learn that network correlation causes the threshold to decrease with the increase of µ, while dynamic correlation has the opposite effect.
We know that the average degree can regulate the strength of the dynamic correlation.In Fig. 2(a), we significantly reduce the average degree, at which point the correlation of the dynamic is greatly enhanced and dominant, so we see the threshold monotonically increases as µ increases.In contrast, Figure 2(b) shows that the threshold decreases monotonically as µ increases, as the average degree increases substantially.In addition, we can see in Fig. 2(b) that the results from Eq. ( 10) match well with those obtained from Monte Carlo simulations.
As a final analysis, we plot the epidemic prevalence ρ I as a function of β/µ for different µ in Fig. 3.For the analytic results in Fig. 3(b), we can see that ρ I appears as a crossing behavior as a function of β/µ, which are the results of considering only the network correlation and ignoring the dynamic correlation.The crossing point occurs precisely at ρ I = 0.5 and the corresponding β/µ = [m⟨a/(a + b)⟩] −1 .For the simulation results of Fig. 3(a), we see that different curves no longer intersect precisely at one point, because the effects of dynamic correlation are not the same for different µ.Comparing the curves with high µ in Fig. 3(a) and Fig. 3(b), such as µ = 10, we can see that ρ I has a large drop for any β/µ, which means that the correlation of dynamic play a key role for high µ.In addition, we can see in Fig. 3 that our Eq.( 8) can be used to predict the prevalence for low µ.
In summary, we investigated the effects of network correlation and dynamic correlation on SIS epidemics of activitydriven networks.In particular, we analyze the effects of the two types of correlations in isolation by comparing theoretical and simulation results.As recovery rate µ increases, the threshold (β/µ) c decreases for the effect of network correlation but increases for the effect of dynamic correlation.Because of this competitive relationship, the curves of the threshold show three types of behavior, which are monotonically increasing, monotonically decreasing, and first increasing and then decreasing but increasing again.Meanwhile, the effect of network correlation means that the prevalence produces a crossing behavior, and the effect of dynamic correlation substantially reduces the prevalence.Finally, our theory can predict the threshold when the average degree ⟨k⟩ is high and predict the prevalence for low recovery rate µ.As a final remark, we hope our study will open the door for other analytical studies of the full, two-dimensional parameter space of the SIS model on temporal networks.ual becoming inactive is λ i (τ) = b.Summing up all of them, we yield the total transition rate ω(τ) = i λ i (τ).
(iii) Time is incremented by dτ = 1/ω(τ).The individual whose state is chosen to change at time τ + dτ is sampled with a probability proportional to λ i (τ).The selected individual changes its state.When the selected individual is activated, it randomly selects m individuals to generate links.When the selected individual becomes inactive, it deletes all links to itself.
(iv) Repeat steps (ii)-(iii) until the number of active individuals and the average degree of the network are stable.
(v) Set t = 0, I 0 initial infected individuals are randomly selected.
(   2), Eq. ( 4), and Eq. ( 5), respectively.Parameters: In addition, we present some simulation results in Fig. S5.Note that the simulation results include the correlation of dynamic.As can be seen from Fig. S5(a), for most nodes (unshaded area), there is little difference in ρ UI a near the threshold even if the dynamic correlation is included.
Finally, by inserting Eq. (S4) into Eq.(S3a), we can obtain Eq. ( 8): We define ρ I a (ρ S a ) as the fraction of individuals in infected (susceptible) state and class a.According to Eq. (S1a), we have dρ I a dt = − µρ I a + βm ρ US a ρ AI a ′ + 2ρ AS a ρ AI a ′ + ρ AS a ρ UI a ′ da ′ , (S1a) = − µρ I a + βm ρ S a ρ AI a ′ + ρ AS a ρ I a ′ da ′ .(S1b) We define ρ I and ρ S as the number of individuals in the infected state and in the susceptible state, i.e., ρ I = ρ I a da and ρ S = ρ S a da.By integrating both sides of Eq. (S1a) with respect to a, we get
FIG.S5.The ρ UIa as a function of a in the steady state near and far from the epidemic threshold in (a) and (b), respectively.Parameters:N = 10 5 , b = 1, m = 3, η = 28, F(x) ∝ x −2.1 with x ∈ [10 −3 , 1].Each simulation runs to t = 5000.The results represent the average of the last 1000 time intervals over 50 independent simulations.The total fraction of nodes in the shaded area is less than one percent.
− b − µ ρ AI − 2βm ρ AI 2 + βm a a + b + ⟨a⟩ ρ UI − βmρ UI ρ AI , + b ρ AI − (⟨a⟩ + µ) ρ UI − βmρ UI ρ AI .(S5b)Appendix D: A theoretical solution for the case of 50% prevalence dρ I dt = −µρ I + βm ρ S ρ AI + ρ AS ρ I .(S2)Consider the special steady state ρ I = ρ S = 0.5, we have vi) At any time t, we calculate each individual's transition rates of the network state λ i 1 (t) and those of the dynamic state λ i 2 (t).The rate for any inactive individual becoming active is λ i 1 (t) = a i .The rate for any active individual becoming inactive is λ i 1 (t) = b.The rate for any susceptible individual becoming infected is λ i 2 (t) = βk inf and k inf is the number of infected neighbors of the focal individual.The rate for any infected individual recovering is λ i 2 (t) = µ.Summing up all of them, we yield the total transition rate ω(t) = i [λ i 1 (t) + λ i 2 (t)].(vii)Time is incremented by dt = 1/ω(t).The individualwhose state is chosen to change at time t + dt is sam-When the network state of the selected individual changes, the rules for generating or disconnecting links are the same as those in step (iii).