key: cord-0010881-90379pjd
authors: Saeedian, M.; Khalighi, M.; Azimi-Tafreshi, N.; Jafari, G. R.; Ausloos, M.
title: Memory effects on epidemic evolution: The susceptible-infected-recovered epidemic model
date: 2017-02-21
journal: Phys Rev E
DOI: 10.1103/physreve.95.022409
sha: 345410ddcf907a621352c077e2494ebfe2a6a367
doc_id: 10881
cord_uid: 90379pjd

Memory has a great impact on the evolution of every process related to human societies. Among them, the evolution of an epidemic is directly related to the individuals' experiences. Indeed, any real epidemic process is clearly sustained by a non-Markovian dynamics: memory effects play an essential role in the spreading of diseases. Including memory effects in the susceptible-infected-recovered (SIR) epidemic model seems very appropriate for such an investigation. Thus, the memory prone SIR model dynamics is investigated using fractional derivatives. The decay of long-range memory, taken as a power-law function, is directly controlled by the order of the fractional derivatives in the corresponding nonlinear fractional differential evolution equations. Here we assume “fully mixed” approximation and show that the epidemic threshold is shifted to higher values than those for the memoryless system, depending on this memory “length” decay exponent. We also consider the SIR model on structured networks and study the effect of topology on threshold points in a non-Markovian dynamics. Furthermore, the lack of access to the precise information about the initial conditions or the past events plays a very relevant role in the correct estimation or prediction of the epidemic evolution. Such a “constraint” is analyzed and discussed.

The study of epidemiology, concerning the dynamical evolution of diseases within a population, has attracted much interest during recent years [1] . Mathematical models of infectious diseases have been developed in order to integrate realistic aspects of disease spreading [2] [3] [4] [5] . A simple and commonly studied model, introduced by Kermack and McKendrick, is the susceptible-infected-recovered (SIR) model [6] . In this model, populations can be in each of three states: susceptible, infected, and recovered (removed), denoted by S, I, and R, respectively. Originally, it is assumed that susceptible individuals become infected with a rate proportional to the fraction of infected individuals in the overall population (fully mixed approximation) and infected individuals recover at a constant rate. The epidemic process presents a (percolation) transition between a phase, in which the disease outbreak reaches a finite fraction of the population, and a phase with only a limited number of infected individuals [7, 8] . The model has also been investigated for population on lattices (e.g., [9] [10] [11] ) or on networks (e.g., [12] [13] [14] ).

For simplicity, we will keep the "medical epidemic" vocabulary hereafter. However, the model has also been interesting for describing nonmedical epidemics, such as for financial bubbles [15, 16] , migration [17] , opinion formation [18, 19] , or internet "worm propagation" [20] [21] [22] . SIR models with distributed delay and with discrete delay have also been studied [23, 24] .

In the usual SIR model, it is assumed that all contacts transmit the disease with the same probability. Moreover, the transmission and recovery coefficients are constant. Hence the state of system at each time does not depend on the previous history of the system: it is a memoryless, so-called Markovian, process. However, real surveys show evidence of a non-Markovian spreading process [25, 26] in agreement with common expectation. The epidemic processes evolution and control, in human societies, cannot be considered without any memory effect. When a disease spreads within a human population, the experience or knowledge of individuals about that disease should affect their response [5] . If people know about the history of a certain disease in the area where they live, they use different precautions, such as vaccination, when possible. Thus, some endogenous controlled suppression of the spreading is expected, although other factors can help [27] [28] [29] . However, knowledge about the history of a disease does not have the same influence at all times. Experience about the prevalence of a disease and precautions related to the "old times" are not always applicable or recommended, hence people tend to follow new strategies against the diseases. In other words, memory of the earlier times could have less effect on the present situation, as compared to more recent times. It can be expected that long-range memory effects decay in time more slowly than an exponential decay, but can typically behave like a power-law damping function.

While much effort has been made so far to determine exact epidemic thresholds in Markovian epidemic models [30] [31] [32] [33] [34] , few works have been devoted to study the non-Markovian aspects of epidemic processes [35, 36] . Furthermore, in this work we focus on long-range memory effects, which means arbitrarily long history can be included. That is in contrast to short-term memory effects, which have been extensively studied. For instance, Dodds and Watts [37] introduced a general model of contagion considering memory of past exposures to a contagious influence. The authors have argued that their model can fall into one of three universal classes, due to the behavior of fixed point curves. Also, in [38] [39] [40] , the authors consider "implicit memory" by applying asynchronous adapting in disease propagation. They show that this type of memory can lead to a first-order phase transition in outbreaks, thus hysteresis can arise in such models [40] .

It is here briefly recalled that fractional calculus is a valuable tool to observe the influence of memory effects on the dynamics of systems [41] [42] [43] [44] [45] [46] , and has been recently used in epidemiological models [47] [48] [49] [50] . Typically, the evolution of epidemiological models is described with differential equations, the derivatives being of integer order. By replacing the ordinary time derivative by a fractional derivative, a time correlation function or memory kernel appears, thereby making the state of the system dependent on all past states. Thus, it seems that such a method based on derivatives with noninteger order, as introduced by Caputo for geophysics problems [51] , is a very proper formalism for such non-Markovian problems. Moreover, Caputo's formalism provides the advantage that it is not necessary to define the fractional order initial conditions, when solving such differential equations [51] [52] [53] [54] . Furthermore, the time correlation function, in the definition of Caputo fractional derivative, is a power-law function, which is flexible enough to reflect the fact that the contribution of more early states is noticeably less relevant than the contribution of more recent ones on the present state of the dynamical system.

Most of the previous works have studied the epidemiological models with fractional order differential equations from a mathematical point of view. They mainly focused on presenting effective mathematical methods in order to solve the corresponding differential equations [49, [55] [56] [57] . For instance, in [57] a mathematical tool (the multistep generalized differential transform method) is introduced to approximate the numerical solution of the SIR model with fractional differential equations. Also in [48] the authors use fractional order differential equations for epidemic models and concentrate on the equilibrium points of the models and their asymptotic stability of differential equations of fractional order. Other variations of the SIR model with fractional derivatives have also been studied. For instance, Seo et al. introduced the SIR epidemic model with square root interaction of the susceptible and infected individuals and discussed the local stability analysis of the model [56] . Also in [49] , numerical solution of the SIR epidemic model of fractional order with two levels of infection for the transmission of viruses in a computer network has been presented.

In all previous works, the authors rarely discuss the effect of fractional order differential equations and memory on the epidemic thresholds and the macroscopic behavior of epidemic outbreaks. Hence, one question remains; we address it in this paper: How does the system robustness change if memory is included in the SIR model? We also use the fractional differential equations, describing the SIR model on structured networks, to see the effect of topology on the evolution of the SIR model including memory effects.

Furthermore the lack of access to accurate information on initial conditions sometimes leads to doubt about epidemic evolution predictions [58] . The same type of difficulty occurs in related problems, such as in opinion formation [59, 60] . Moreover, it may also happen in certain cases that individuals do not believe in old strategies in order to avoid the disease.

This means that the initial time for taking into account the disease control memory is shifted toward more recent times: thereafter, the dynamics is evolving with a new fraction of susceptible and infected individuals, different from that predicted by the solution of the differential equations. In contrast, the fractional calculus method allows us to choose any arbitrary initial time at which the effect of initial conditions can be introduced on the spreading dynamics with a memory content. The interest of fractional calculus will appear through such aspects in the core of the paper.

Thus, the paper is organized as follows. In Sec. II, following Caputo's approach, we convert the differential equations of the standard SIR model to the fractional derivatives, thereby allowing us to consider memory effects. Using numerical analysis results (Sec. III), we discuss the influence of memory on the epidemic thresholds in Sec. III A. We also discuss the dynamics of a non-Markovian epidemic process, when choosing different initial conditions or modifying the proportions of agents at a given time in Sec. III B. To complete our discussion, we study the dynamics of the model on structured networks in Sec. IV. We also point out that we have observed qualitatively similar results for the SIS (susceptible-infectedsusceptible) epidemic model. The conclusions are found in Sec. V.

The evolution of the standard SIR model is described by a set of coupled ordinary differential equations for susceptible (S), infected (I ), and recovered(R) individuals, respectively given by

in which β and γ are infection and recovery coefficients, respectively. The infected individual makes β contacts per unit time producing new infections within a mean infectious time of order 1/γ . The evolution of the model is controlled by quantity β/γ , such that above the epidemic threshold, (β/γ ) c , the disease spreads among a finite fraction of individuals. These (ordinary) differential equations describe a Markov epidemic process, in which the state of individuals at each time step does not depend on previous steps. The set of Eqs. (1) can be solved iteratively until time t. In particular, the fraction of susceptible individuals at time t, denoted as S t , can be determined. In fact, 1 − S t is the size of outbreaks, i.e., the population that has or has had the disease until time t.

In order to observe the influence of memory effects, first we rewrite the differential equations (1) in terms of time-dependent integrals as follows:

in which κ(t − t ) plays the role of a time-dependent kernel and is equal to a delta function δ(t − t ) in a classical Markov process. In fact, any arbitrary function can be replaced by a sum of delta functions, thereby leading to a given type of time correlations. A proper choice, in order to include long-term memory effects, can be a power-law function which exhibits a slow decay such that the state of the system at quite early times also contributes to the evolution of the system. This type of kernel guarantees the existence of scaling features as it is often intrinsic in most natural phenomena. Thus, let us consider the following power-law correlation function for κ(t − t ):

in which 0 < α 1 and (x) denotes the gamma function. The choice of the coefficient 1/ (α − 1) and exponent (α − 2) allows us to rewrite Eqs. (2) to the form of fractional differential equations with the Caputo-type derivative. If this kernel is substituted into Eqs. (2), the right-hand side of the equations, by definition, are fractional integrals of order

Applying a fractional Caputo derivative of order α − 1 on both sides of each Eq. (2), and using the fact the Caputo fractional derivative and fractional integral are inverse operators, the following fractional differential equations can be obtained for the SIR model:

t denotes the Caputo derivative of order α, defined for an arbitrary function y(t) as follows [51] :

Hence, the fractional derivatives, when introducing a convolution integral with a power-law memory kernel, are useful to describe memory effects in dynamical systems. The decaying rate of the memory kernel (a time-correlation function) depends on α. A lower value of α corresponds to more slowly decaying time-correlation functions (long memory). Hence, in some sense, the strength (through the "length") of the memory is controlled by α. As α → 1, the influence of memory decreases: the system tends toward a memoryless system. Note that for simplicity, we assume the same memory contributions (same value of α) for different states of S, I , and R. Obviously, more complicated functions than Eq. (3) and taking into account different α i (i = 1,2,3) could be investigated in further work to take into account different time scales.

Although analytical solutions of Eqs. (4) are hard to obtain for the general case, they can be obtained at the early stage of the epidemic under a linearization approximation. In this case, it turns out that the number of infected individuals behaves as a Mittag-Leffler function [52] :

in which ζ is a constant, which depends on the initial conditions [52] . In particular, for α = ζ = 1, the Mittag-Leffler function is the exponential function. Thus, in the early stage of epidemic dynamics, the growth rate of the infected population in Eq. (6) is positive, if β − γ > 0. Therefore, the number of infected individuals grows exponentially in such a case, for β > γ , as of course is expected for the standard memoryless SIR model. The same reasoning applies in order to determine the epidemic threshold for α < 1.

Let it be reemphasized that Eqs. (4) consist in a system of coupled nonlinear differential equations of fractional order, in the following general form:

where i = 1,2,3 and y (1) ,y (2) ,y (3) denote S,I,R cases, respectively. Also, y (i0) are constants which indicate the initial conditions.

To solve the equations, we use the predictor-corrector algorithm, which is well known for obtaining a numerical solution of first-order problems [61] [62] [63] . It is assumed that there exits a unique solution for each of y (i) on the interval [0,T ] for a given set of initial conditions. Considering a uniform grid {t n = nh : n = 0,1,2, . . . ,N}, in which N is an integer and h ≡ T /N, each Eq. (7) can be rewritten in a discrete form,

where the coefficients b n−k−1 refer to the contribution of each of the n − 1 past states on the present state of n. The coefficients are given by

Thereby after solving Eq. (7), numerically, the influence of memory on the evolution of the SIR epidemic model can be analyzed. As mentioned in the Introduction, let us consider two pertinent aspects successively: the finite time behavior and the role of changing initial conditions.

Let us compare the evolution of a system including memory effects with the memoryless case. We solve Eq. (7) with initial conditions y (10) = S 0 = 1 − , y (20) = I 0 = . Figure 1 is zero (with accuracy 10 −4 ) for small values of β/γ . The specific value of β/γ , in which the epidemic size starts to get a nonzero value, is identified as the epidemic threshold point.

The stationary time for a memoryless system (α = 1) is t = 100. With decreasing the value of α (including memory) the system needs much time to reach the stationary state. Hence at t = 100, the threshold point is shifted to the higher value of β/γ . Figure 2 shows that the threshold point is increased with decreasing of α for a finite time t. Furthermore, as can be seen in Fig. 1 , the size of outbreaks decreases for decreasing α.

Let the interval [t 0 ,t] be the time interval in which memory effects are taken into account. In Fig. 3 , we compare the evolution of the model with memory for different values of the finite time t. The memory effects are considered for a weight α = 0.2. It is seen that as time evolves the influence of memory decreases, since memory effects decay in time like Fig. 4 . Furthermore, for a given β/γ value, it appears that there is more time available for disease spreading, whence more individuals become infected.

Recall that the dynamics of a non-Markovian process is directly influenced by all events from the beginning of the process. However, some loss of information about some period of time in the past may lead one to consider that the influence of memory might not need to be considered as continuous. It may happen, in many social networks, that individuals do not have enough information about the history of a disease, as recent cases and studies indicate; e.g., see [64] [65] [66] [67] . Only after several individuals have already been infected, do people start to increase their knowledge about the disease and take different precautions. The question arises on how the "initial time" at which a non-Markovian process is started affects the subsequent dynamics of the process.

If two Markovian processes start at two different times, the evolution of both processes is identical. However, the scenario is quite different for a non-Markovian process, i.e., in which the memory plays a role. This is illustrated through Fig. 5 where the fractions of susceptible, infected, and removed individuals are compared in the case of two Markov and non-Markov epidemic processes. Continuous and dashed (black) lines correspond to a system with and without including memory effects, respectively, evolving from the same initial time t = 0. As can be seen, the fraction of susceptible individuals is greater in a system with inclusion of memory effects with respect to that ignoring the memory [ Fig. 5(a) ]. In other words, the experience and knowledge which individuals have about the disease are obviously helping them to protect themselves against the disease. Equivalently, in a system including the memory effects, the infection grows more slowly as seen in Fig. 5(b) .

Thereafter, consider that a non-Markovian process, including memory effects, has evolved until a specific time t 1 . Let the process be continuing its evolution, but let the memory of the system be removed at that time. This corresponds to having a new initial time and new initial conditions for the epidemics spreading. The process can be continued without or with memory. The Markovian case is trivial thereafter and thus not discussed. Instead, consider that memory effects are only taken into account at this starting "new initial time." In other words, let the population ignore the disease control history (memory) until t 1 ; let the system continue its dynamics but taking into account memory effects thereafter from t 1 . The initial conditions for the evolution of the system are now a fraction of susceptible and infected individuals at time t 1 . The curves with square symbols in Fig. 5 correspond to what happens for different "new initial times" t 1 = 30,70, for the dynamics of such a non-Markovian epidemic process. As can be seen, at the beginning of the dynamics, the fraction of susceptible individuals is reduced, since people do not know about the disease. However, as soon as it is influenced by memory, the system becomes more resilient to the spreading. Hence, the fraction of S individuals remains greater as compared to that with a memoryless system, having started at t = 0. In a similar manner, the fraction of infected and removed individuals deviate from the original one and tend toward the populated states of a memoryless system when the memory from further past times is included. In this case, the curves become closer to the dashed curve corresponding to a memoryless system. That means that the system loses the information related to past times and tends to present a behavior similar to a memoryless system. Finally, one can consider "to remove the memory" of an epidemic process at various times. At each time step, the system is supposed to lose (or practically negate) the information about the disease before some "reawareness time" (see also [68] ) and to continue its dynamics regardless of the past. For illustration, consider the case of such a sudden Notice that in this particular illustrative case, the behavior of the system is seen to be close to the dynamics of a memoryless system, since contributions of the memory of the system are sometimes removed. Such an illustration points to the interest of the model in order to compare it with the case of epidemics spreading waves [65] ; for completeness, let it be pointed out that the connection of periodic epidemics to SIR models has been already mentioned [69] : flu is yearly recurrent. Notice also that the value of α could be modified at each new awareness time, but this investigation goes outside the present paper.

So far we have considered the fully mixed approximation, such that an infected individual is equally likely to spread the disease to any other individual. However, in the real world an individual connects to a small fraction of people. Hence, as is well known, more realistic modeling can be studied through networks, where their topology has a significant effect on the epidemic process [70] [71] [72] . For homogeneous networks, each individual has the same number of connections k ≈ k and disease propagates with spreading rate β k . In this case, it is obvious that the epidemic threshold ( β γ ) c is simply replaced by ( β γ ) c k . It is also true for the case of fractional differential equations (4) . In other words, threshold point in Fig. 1 for each value of α is shifted to ( β γ ) c k . Now, let us consider heterogeneous scale-free networks with degree distribution P (k) ∼ k −λ . In heterogeneous meanfield approximation, it is assumed that all nodes are statistically equivalent and thus one can consider groups of nodes with the same degree k. With this assumption the ordinary differential equations describing the SIR model are given by where

and i k , s k , and r k denote the density of infected, susceptible, and removed nodes in each group, respectively. It turned out that in scale-free networks characterized by a degree exponent 2 < λ 3, there is no epidemic threshold [72] . Following the same procedure presented in Sec. II, we can rewrite Eqs. (10) to the fractional derivatives, as follows:

, For a network with degree exponent λ = 3, we solve Eqs. (12) numerically. Figure 6 shows the evolution of the fraction of total infected individuals i(t) = k P (k)i k (t), with considering memory effects with different values of α. While for a memoryless SIR model (α = 1), the system reaches a stationary state after a short time (t 20), the stationary time is increased with decreasing the value of α. Furthermore, we obtain the size of outbreaks at a finite time. Figure 7 shows 1 − S t , measured with accuracy 10 −5 until t = 100 for different values of α. As we can see, the epidemic threshold is always zero, as it is for Markov epidemic spreading on scale-free networks with λ = 3. However, the size of epidemic decreases with decreasing α. The same results are obtained for networks with 2 < λ < 3. However, for λ > 3, the epidemic threshold is shifted with including the memory, similar to what is observed for the homogeneous networks.

Memory plays a significant role in the evolution of many real dynamical processes, including the cases of epidemic spreading. Here we have reported a study on the evolution of the SIR epidemic model, considering memory effects. Using the fractional calculus technique, we show that the dynamics of such a system depends on the strength of memory effects, controlled by the order of fractional derivatives α. At finite times, including memory effects, the epidemic threshold (β/γ ) c is shifted to higher values than those for memoryless systems, at values depending on the memory decay rate α. In the case that the model evolves on heterogeneous scale-free networks with 2 < λ 3, the threshold point is always zero. However, the fraction of individuals who are infected or recovered, is reduced if the memory "length" increases. Hence, memory renders the system more robust against the disease spreading. If the epidemic process evolves further in time, for a fixed memory strength, (i) the disease can infect more individuals and (ii) the epidemic threshold is shifted to smaller values and tends to the memoryless case values.

Furthermore, we have shown the following result: the evolution of an epidemic process, including memory effects, much depends on the fraction of infected individuals at the beginning of the memory effect insertion in the evolution. During a non-Markovian epidemic process, if the system abruptly loses its memory at a definite time and if from that time on, one lets the non-Markovian process continue again, starting with the number of infected individuals at that time, the dynamics of the system deviates from the basic case, in which the system continuously includes memory effects from the beginning of the process.

Our observations are obtained from a simple epidemiological model: the SIR model. Obviously many parameters are here assumed to be constant. We are aware that some, e.g., policy, feedback might influence the parameter values. They may depend on space, groups, and time. External field conditions may also surely influence real aspects. However, we guess that many qualitative behaviors as those presented here are likely to be quite generally found in reality. More advanced epidemic models, based on various types of complex networks, are surely interesting subjects for further investigations, in line with investigations such as, e.g., in [30] [31] [32] [33] [34] . We also wish to point out that we have observed qualitatively similar results for the SIS epidemic model. Finally, we may claim that our results are not limited to the epidemiological ("medical") models but also can be extended for analogous epidemic spreading of rumors, gossip, opinions, religions, and other topics pertinent to epidemics on many social networks.

Knowledge Epidemics and Population Dynamics Models for Describing Idea Diffusion, Models of Science Dynamics: Encounters Between Complexity Theory and Information Sciences

Introduction to Percolation Theory

The Logistic Map and the Route to Chaos: From the Beginning to Modern Applications

A mixed abstraction level simulation model of large-scale internet worm infestations

Measurement and Analysis of Worm Propagation on Internet Network Topology, Proceedings of the 13th International Conference on Computer Communications and Networks (ICCCN, 2004) (unpublished)

Climate Change and Crop Production, CABI Climate Change Series

Fractional Calculus: An Introduction for Physicists

Applications of Fractional Calculus in Physics

Fractional Differential Equations

Theory and Applications of Fractional Differential Equations, North-Holland Mathematics Studies

The FracPECE Subroutine for the Numerical Solution of Differential Equations of Fractional Order, Forschung und wissenschaftliches Rechnen

Changing Waves: The Epidemics of 1832 and 1854

Proc. Natl. Acad. Sci. USA

It could be instructive to study fractional order operators within a geometric interpretation (see interesting references in [53] ). Here we compare the time scales in fractional-and integer-order dynamics. To image a geometric interpretation, let us consider the right-sided fractional integral of order α,and write it in the formwhereIf we compare Eq. (A2) with its counterpart x α t = t 0 v(τ )dτ , in the homogeneous time scheme, the main difference is related to the different time variables T and τ . Notice that time variable T t (τ ) has a scaling property. If we take t 1 = kt and τ 1 = kτ , then T t 1 (τ 1 ) = k α T t (τ ). Hence, in the fractional order dynamics, the time is "accelerating" in the early time and after that it is "slowing down," as sketched in Fig. 8 .In this case, the "passing time" in the two axes of time is not the same. For this reason, in epidemic "fractional" dynamics, the epidemic threshold is shifted to the higher values. A lower α indicates a "stronger" (long-lasting) memory and a more pronounced shift of threshold point. However, if one waits long enough, the same behavior will be observed in both fractional and usual homogeneous time. In fractional dynamics, after a "very long" time, the threshold point coincides with the one appearing in integer-order dynamics.