There has been growing interest in exploring the interplay between epidemic spreading with human response, since it is natural for people to take various measures when they become aware of epidemics. As a proper way to describe the multiple connections among people in reality, multiplex network, a set of nodes interacting through multiple sets of edges, has attracted much attention. In this paper, to explore the coupled dynamical processes, a multiplex network with two layers is built. Specifically, the information spreading layer is a time varying network generated by the activity driven model, while the contagion layer is a static network. We extend the microscopic Markov chain approach to derive the epidemic threshold of the model. Compared with extensive Monte Carlo simulations, the method shows high accuracy for the prediction of the epidemic threshold. Besides, taking different spreading models of awareness into consideration, we explored the interplay between epidemic spreading with awareness spreading. The results show that the awareness spreading can not only enhance the epidemic threshold but also reduce the prevalence of epidemics. When the spreading of awareness is defined as susceptible-infected-susceptible model, there exists a critical value where the dynamical process on the awareness layer can control the onset of epidemics; while if it is a threshold model, the epidemic threshold emerges an abrupt transition with the local awareness ratio α approximating 0.5. Moreover, we also find that temporal changes in the topology hinder the spread of awareness which directly affect the epidemic threshold, especially when the awareness layer is threshold model. Given that the threshold model is a widely used model for social contagion, this is an important and meaningful result. Our results could also lead to interesting future research about the different time-scales of structural changes in multiplex networks.

Network modeling plays an important role in identifying statistical regularities and structural principles common to many complex systems. As a representative way to describe the rapidly changing interactions among individuals, time varying network is gaining more and more interest. Hence, in order to explore the relationship between two processes accounting for the spread of epidemics, and the spread of awareness information to prevent its infection, we propose a multiplex network model with two layers, of which the information spreading layer is a time varying network and the contagion layer is a static network. By extending the microscopic Markov chain approach (MMCA), we are able to derive the epidemic threshold analytically. The numerical results show high accuracy compared with extensive Monte Carlo (MC) simulations, regardless of whether the information spreading process is defined as the susceptible-infected-susceptible model (SIS) model or the threshold model. Specifically, when the spreading of awareness is defined as the SIS model, there exists a critical value where the dynamical process on the awareness layer can control the onset of epidemics; while if it is a threshold model, the epidemic threshold emerges an abrupt transition with the local awareness ratio α approximating 0.5. The results also show that changes in the topology of the awareness layer could change the location of the meta-critical point, which is particularly important since a slight delay or advance of the epidemic could have important implications to vaccination strategies. In sum, our results give us a better understanding of the coupled dynamical processes and highlight the importance of taking the spread of information into consideration in the control of epidemics.

The problem of modeling how epidemics spread among individuals has been extensively studied by people from various fields for many years.1–5 In the field of physics, most approaches to these problems are related to the theory of phase transition, statistical physics, and critical phenomenon.6–10 Especially in last decades, the study of complex networks has provided new grounds to the understanding of epidemic dynamics for physicists.11,12 They proposed various models, such as the classical susceptible-infected-susceptible model (SIS),13 susceptible-infected-recovery model (SIR),14 and so on,15,16 to shed light on modeling the epidemic dynamics. Different models have addressed the problem from different perspectives, e.g., the frequency of contacts between people,17 duration of the disease,18 and risk perception.19 

Furthermore, in network modeling of a single layer topology structure, the standard approach assumes that all links have the same activation time, which is clearly an abstraction of any real topology structure and cannot reveal the true relationships in the complex system very well.20 Therefore, as a natural way to describe the interrelated complex connections among individuals, the multiplex network is gaining more and more attention. It can overcome the drawbacks of the simplex network in which individuals interact only through one network, while in reality (same as in epidemic and information spreading), the same set of nodes might interact through multiple networks with different topologies and dynamical properties.21–25 Each layer could have particular features that are different from the rest, which make the multiplex network define a richer structure of interactions.26 Under the framework of the multiplex network, there has been growing interest in exploring the interplay between epidemic spreading with human response,27–33 since it is natural for people to take measures when they become aware of epidemics. Two interdependent layers can be used for modeling different dynamical processes, and the interactions between layers represent the coevolution of the processes. Several works have studied the problem by changing the model of information spreading, such as by considering the information spreading process as a SIS model; Granell et al.27 discovered the emergence of a metacritical point where the diffusion of awareness is able to control the onset of epidemics, whereas Guo et al.28 replaced the SIS model by a threshold model and found a two-stage effect on epidemic threshold. Besides, Granell et al.29 also analyzed the consequences of a massive broadcast of awareness (mass media) on the final outcome of the epidemic incidence, and the results showed that mass media could make the metacritical point disappear.

Until now, people pay much attention on the spreading model for the multiplex network, while the topology structure of it has gained little attention. In other words, all these models can be defined as connectivity model, where connections among individuals are long-lived elements.34–37 However, with the development of various kinds of social networks, in many cases, the interactions among individuals are rapidly changing and are characterized by processes with duration being defined on a very short time scale, such as information diffusion over an email, mobile calls among individuals, and so on.38 Since a time-aggregated representation of network's interactions neglects the time-varying nature of real systems,35 it is meaningful for us to model the information layer as a time varying network and to study the effects that information spreading process has on epidemic spreading. This scenario is common in reality if we consider the spreading of epidemic information originating from social media, for example, the spreading of Severe Acute Respiratory Syndromes. During this period, people can know about this epidemic from many sources, such as television (TV), radios, and Internet. Let us take TV as an example: the announcer can be regarded as an activated node and the connections between him and listeners are time-varying, which makes the proper topology of the connections be a time-varying network. At the same time, in order to decrease the infection probability, everyone needs to avoid contacts with strangers. And then, it is obvious that the contagion spreading layer can be treated as a static one. Therefore, in this paper, to fill this gap, a multiplex network with two layers is built. Specifically, the information spreading layer is a time varying network generated by the activity driven model,34 and the other one is a static network. We extend the microscopic Markov chain approach (MMCA)27 to derive the epidemic threshold of the model analytically. Compared with extensive Monte Carlo (MC) simulations, the method shows high accuracy for the prediction of the epidemic threshold. Besides, by considering different spreading models, namely, the SIS model and the threshold model, of the time varying network, we explored the interplay between epidemic spreading with information spreading. The results show that the coexistence of antagonistic spreading effects raises interesting physical phenomena, such as not only makes the epidemic threshold become larger but also lowers the density of the final epidemic size. Corresponding to different models, the metacritical point and two-stage effects also exist in the time varying cases, which verify again that these phenomena are rooted in the competition principle of the multiplex network model. Furthermore, we also find that the temporal changes in the topology hinder the spread of awareness which directly affect the epidemic threshold, especially when the awareness layer is the threshold model.

Let us start by defining the specific setup of the multiplex network we analyze. The awareness layer corresponds to the spread of information about epidemics, where individuals can be aware or unaware. Besides, there is a certain probability for the two states transforming between each other. More specifically, aware individuals have a probability of δ forgetting epidemics while unaware individuals can become aware with the probability λ through communication with aware individuals. As an important way to describe information diffusion process, the activity driven model is used to generate the connections among individuals on the awareness layer.34 Meanwhile, on the epidemic layer, there are two categories of individuals, i.e., S (susceptible) individuals and I (infected) individuals. Furthermore, the contagion layer, where the epidemic spreading process described by the classic SIS model occurs, is a static one. Because of the existence of coupled dynamical processes, it produces three kinds of states for each individual, including US (unaware and susceptible), AS (aware and susceptible), and AI (aware and infected), as illustrated in Fig. 1. Under the assumption of the SIS model for the spreading of epidemic, a susceptible individual can be infected by an infected neighbor with the probability β, and at the same time, an infected individual can recover to be susceptible with probability μ. If an individual is infected, it is natural that this individual becomes aware of the epidemic. However, if one susceptible individual is aware of the epidemic, the infectivity β can be reduced by a factor, which results from individual's behavioral responses to the presence of the epidemic.27 We use βU and βA to represent the infection rates without and with awareness, respectively. For the sake of simplicity, we will present the results for βA=0 in this paper, which correspond to complete immunity for aware individuals.

FIG. 1.

A simple example of the structure of the multiplex network proposed in the paper. The awareness layer is a time varying network where the spreading of awareness happens. At each time step, the topology structure of the awareness layer is built according to the activity driven model. Relying on the newly built network, the individuals interact with each other to change the states: unaware (green node) and aware (yellow node). The other layer corresponds to the network where epidemic transfers among two kinds of nodes, namely, susceptible (green node) and infected (red node). Individuals on two layers are the same. In particular, only three kinds of states exist in the multiplex network: unaware and susceptible, aware and infected, and aware and susceptible. Note that in the figure, red individuals on the time varying network represent the activated ones and each activated individual can have 2 connections with other individuals at each time step.

FIG. 1.

A simple example of the structure of the multiplex network proposed in the paper. The awareness layer is a time varying network where the spreading of awareness happens. At each time step, the topology structure of the awareness layer is built according to the activity driven model. Relying on the newly built network, the individuals interact with each other to change the states: unaware (green node) and aware (yellow node). The other layer corresponds to the network where epidemic transfers among two kinds of nodes, namely, susceptible (green node) and infected (red node). Individuals on two layers are the same. In particular, only three kinds of states exist in the multiplex network: unaware and susceptible, aware and infected, and aware and susceptible. Note that in the figure, red individuals on the time varying network represent the activated ones and each activated individual can have 2 connections with other individuals at each time step.

Close modal

In order to analyze the coupled dynamical processes on top of the multiplex network, we need to explore the details of the model. As described above, the multiplex network consists of two layers with different topology structures. The awareness layer is a time-varying network which is generated by the activity driven model, while the contagion layer is a static network. The activity driven model used for the spreading of information can be depicted as follows:34 first, we consider N individuals and assign to each individual i an activity rate ai=ηxi, which is defined as the probability per unit time to create new contacts with other individuals. The activity potentials xi are assigned according to a probability distribution F(x) and are bounded in the interval ξxi1. Then, at each time step, with probability aiΔt, each individual i is activated and generates m connections with other randomly selected individuals. After we build the network, information spreading process, which is the same as in the classic SIS model, will happen on it. At the next time step, we delete all the edges and repeat the above steps. It is obvious that in the model, individuals do not have memory of previous steps. Since many empirical researches show that power-law distribution of F(x) fits the real world networks,34 in the paper, we make F(x) take the form: F(x)xγ. As a result of the coupled dynamical processes, it is worth noting that each individual in this multiplex network can only be in one of the three kinds of states: unaware and susceptible (US), aware and infected (AI), and aware and susceptible (AS). To help readers have a clear understanding of our model, we list the meanings of some key parameters used in our model in Table I.

TABLE I.

The details of some key parameters.

ParameterDescription
λ Transition probability from unaware to aware 
δ Transition probability from aware to unaware 
βU Infection rate for unaware individual 
βA Infection rate for aware individual 
μ Recovery rate for infected individual 
xi The activity potential of individual i 
ai The activity rate of individual i 
η Rescaling factor of the activity potential 
ξ The lower bounder of individual's activity potential 
m The number of links an active individual generates 
α The local threshold ratio 
ParameterDescription
λ Transition probability from unaware to aware 
δ Transition probability from aware to unaware 
βU Infection rate for unaware individual 
βA Infection rate for aware individual 
μ Recovery rate for infected individual 
xi The activity potential of individual i 
ai The activity rate of individual i 
η Rescaling factor of the activity potential 
ξ The lower bounder of individual's activity potential 
m The number of links an active individual generates 
α The local threshold ratio 

In the following, with the help of MMCA method which can be illustrated by the way of probability tree,27–30 we are able to obtain the epidemic threshold of our model. In Fig. 2, we reveal the possible states and their transitions in our model. Here, let {aij(t)}N×N represents the adjacent matrix of the awareness layer at time t. Meanwhile, {bij}N×N is the adjacent matrix of the contagion layer. Since at time t individual i can only be in one of the three states, we denote the probabilities of the three states as piAI(t),piAS(t),piUS(t), respectively. Furthermore, on the awareness layer, we define the probability for unaware individual i staying unaware as ri(t); on the contagion layer, we separately define the probabilities for individual i not being infected by any neighbors as qiA(t),qiU(t) if i is aware and unaware.

FIG. 2.

The probability tree for the transitions of states. The states include AI (aware and infected), US (unaware and susceptible), and AS (aware and susceptible). In the probability tree, μ represents transition probability from infected to susceptible, δ represents transition probability from aware to unaware. Meanwhile, qA and qU represent the transition probability for individual not being infected by neighbors if it is aware or unaware, respectively. r represents probability for unaware individual staying unaware. The coupled dynamical processes take place consecutive as time goes by.

FIG. 2.

The probability tree for the transitions of states. The states include AI (aware and infected), US (unaware and susceptible), and AS (aware and susceptible). In the probability tree, μ represents transition probability from infected to susceptible, δ represents transition probability from aware to unaware. Meanwhile, qA and qU represent the transition probability for individual not being infected by neighbors if it is aware or unaware, respectively. r represents probability for unaware individual staying unaware. The coupled dynamical processes take place consecutive as time goes by.

Close modal

With respect to the above definitions, considering the coupled dynamical processes we have27 

ri(t)=j=1N(1aji(t)pjA(t)λ),qiA(t)=j=1N(1bjipjAI(t)βA),qiU(t)=j=1N(1bjipjAI(t)βU).
(1)

Hence, according to the probability tree which reveals the transitions of three different states, we can derive a discrete-time version of the evolution of our model by means of Markov chain method, which is

piUS(t+1)=piAI(t)δμ+piUS(t)ri(t)qiU(t)+piASδqiU(t),piAS(t+1)=piAI(t)μ(1δ)+piUS[1ri(t)]qiA(t)+piAS(1δ)qiA(t),piAI(t+1)=piAI(t)(1μ)+piUS(t){[1ri(t)][1qiA(t)]+ri(t)[1qiU(t)]}+piAS(t){δ[1qiU(t)]+(1δ)[1qiA(t)]}.
(2)

For the purpose of calculating the epidemic threshold, it is necessary to explore the stationary solution of the system of Eq. (2). When time t, there exists an epidemic threshold βcU for the coupled dynamical processes, which means the epidemic can outbreak only if ββcU. By letting t, the probabilities of three states piUS,piAS,piAI fulfill the condition that limtpiUS(t+1)=piUS(t+1)=piUS,limtpiAS(t+1)=piAS(t+1)=piAS,limtpiAI(t+1)=piAI(t+1)=piAI. Since around the epidemic threshold βc, the infected probability piAI=ϵi1, the probabilities qiA and qiU can be simplified as qiA(1βAΣjbjiϵj),qiU(1βUΣjbjiϵj), respectively. Therefore, inserting these approximations into Eqs. (2) and omitting higher order items, Eqs. (2) is reduced to the following form:

piUS=piUSri+piASδ,piAS=piUS(1ri)+piAS(1δ),μϵi=piUS((1ri)βAΣjbjiϵj+riβUΣjbjiϵj)+piAS(δβUΣjbjiϵj+(1δ)βAΣjbjiϵj).
(3)

Afterwards, by analysing the probability piAI(ϵi), we can get

μϵi=(piASβA+piUSβU)Σjbjiϵj.
(4)

It is clear that piAI+piAS+piUS=1, where piA=piAI+piAS. Noting that piAI=ϵi1, we get piAIpiA and piUS=1piAIpiAS=1piA. Hence

μϵi=βU(1piA)Σjbjiϵj,
(5)

that is, to say, Eq. (6) can be reduced to

Σj[(1piA)bjiμβUtji]ϵj=0,
(6)

where tji are the elements of the identify matrix. Let H be a matrix whose element hij equals (1piA)bji. Then, it is obvious that the Eq. (6) has nontrivial solutions if and only if μβU is the eigenvalue of matrix H. Consequently, the epidemic threshold βc is the one which satisfies μβU=Λmax, where Λmax is the largest eigenvalue of matrix H, and it is easy to get

βcU=μΛmax.
(7)

It is worthy noting that for single layer time varying network, the epidemic threshold can be calculated as34 

λδ1m1a+a2,
(8)

where a=1Nj=1Naj,a2=1Nj=1Naj2. This means that only for λδm(a+a2)=λc epidemics can outbreak on top of a single layer time varying network. If we consider the two spreading processes separately, therefore, when λ<λc, Eq. (7) is reduced to βc=μΛmax(B) as defined in Ref. 27; thus, (λc,βc) also defines a metacritical point for the epidemic spreading in our model.

Here, to show the validity of the approach discussed above, we have performed extensive MC simulations on multiplex network with 5000 nodes on each layer. The activity driven model, which is used to generate the information spreading layer, is configured as follows: the power law distribution of the activity rate xi satisfies F(x)x3, and the edge number m that one active node can have is 8. Besides, the rescaling factor η is set to be 10 and ξ=103. The other layer is a scale free network, of which the degree distribution P(k)k3. After having these settings, it is easy for us to obtain the boundary of the metacritical point (1m(a+a2),1Λmax(B)), which equals (0.2677, 0.1390) in our case. In Fig. 3, we show the comparison between analytic results calculated by Eq. (7) and MC simulations under various conditions. All the simulations start from a fraction ρ0 of randomly chosen infected nodes and ρ0 is fixed to be 0.2. At each time step, all the neighbours of an infected node become infected with the same probability β and the infected node recovers at a rate μ. The same process fits for the spreading of information—what we need to do is just replace the infected probability β and recovery rate μ with the aware probability λ and forgetting rate δ, respectively. Iterate the rules of the coupled dynamical processes with parallel updating until the density of infected nodes ρI is steady. In order to reduce the fluctuation of the density ρI, we make time average that satisfies ρI=1TΣt=t0t=t0+T1ρI(t) and take T = 100 (t0=901).

FIG. 3.

The comparisons of the epidemic threshold between analytic results (MMCA method) and MC simulations as a function of the infected probability β and the aware probability λ. The green dashed line is the MMCA results for the epidemic threshold and the heatmap represents the density of infected nodes ρI obtained by MC simulations under different parameters. From top left to bottom right, the recovery rate μ and the forgetting rate δ are set to be (a) μ = 0.2, δ = 0.8, (b) μ = 0.4, δ = 0.6, (c) μ = 0.6, δ = 0.4, and (d) μ = 0.8, δ = 0.2, respectively. The dashed rectangle corresponds to the area where the metacritical points are located, which are bounded by the topological characteristics of each layer 1m(a+a2) and 1Λmax(B). The four phase diagrams are obtained by averaging 100 MC simulations for each point in the grid 50 × 50.

FIG. 3.

The comparisons of the epidemic threshold between analytic results (MMCA method) and MC simulations as a function of the infected probability β and the aware probability λ. The green dashed line is the MMCA results for the epidemic threshold and the heatmap represents the density of infected nodes ρI obtained by MC simulations under different parameters. From top left to bottom right, the recovery rate μ and the forgetting rate δ are set to be (a) μ = 0.2, δ = 0.8, (b) μ = 0.4, δ = 0.6, (c) μ = 0.6, δ = 0.4, and (d) μ = 0.8, δ = 0.2, respectively. The dashed rectangle corresponds to the area where the metacritical points are located, which are bounded by the topological characteristics of each layer 1m(a+a2) and 1Λmax(B). The four phase diagrams are obtained by averaging 100 MC simulations for each point in the grid 50 × 50.

Close modal

From the results above, it is clear that the MMCA method has a high accuracy to predict the epidemic threshold, no matter what values other parameters are set to be. Besides, there are also some interesting phenomena revealed by our model. When the value of μ is large, while the value of δ is small, the value of aware probability λ has an obvious effect on the epidemic threshold. The reason is that at this condition, it is easy for the aware nodes to stay in the aware state, while at the same time, the recovery ability of the nodes is strong. These two effects lead the epidemics to be difficult to outbreak. Hence, if we increase λ, which means that the unaware nodes can become aware more easily, the epidemic threshold becomes larger. However, if μ is small and δ is large, it is easy for the epidemics to outbreak. And then, even if we increase λ and nodes become aware with higher probability, the large value of δ makes nodes forget the epidemics in a short time, which will decrease the spread speed of awareness. That is why in this case, λ has a little effect on the epidemic threshold. Furthermore, note that there also exists a region bounded by [1m(a+a2),1Λmax(B)] where the metacritical point is localized, as illustrated in Ref. 27.

In order to explore the difference of epidemic spreading between our model with other models, including the classic SIS model on a single layer network and a static multiplex network, in Fig. 4, we compare these models under different initial conditions. The multiplex network, of which the topology structure of epidemic spreading layer is the same as that of the single layer network, is the one defined in Fig. 3. As shown in Fig. 4, with the help of awareness spreading, not only the epidemic threshold, but also the final density of the infected nodes ρI is smaller than that of the single layer. However, compared with the static multiplex network, it is obvious that the time variability of the awareness layer weakens the suppression effects on the spreading of epidemics, since the epidemic threshold is smaller and the final epidemic size is larger than that of the static multiplex network.

FIG. 4.

MC simulations for epidemic spreading on multiplex network (red lines) and single layer network (blue lines), as well as a two-layer static multiplex network (green lines). The topology structure of the single layer network is the same as the contagion layer of the multiplex network defined above. As for the static network, except that the awareness layer is a static network, its topology structure is totally the same as that of the time-varying network. On each panel, according to different values of the aware probability λ and the forgetting rate δ, we plot four dashed lines for the coupled dynamical processes. The red circle and green diamond line are under the condition when δ = 0.2 and λ = 0.8, while for the red square and green triangle line, δ = 0.8 and λ = 0.2. As for the recovery rate μ, the value is set to be (a) μ = 0.2, (b) μ = 0.4, (c) μ = 0.6, and (d) μ = 0.8. Each line is obtained by averaging 50 independent MC simulations.

FIG. 4.

MC simulations for epidemic spreading on multiplex network (red lines) and single layer network (blue lines), as well as a two-layer static multiplex network (green lines). The topology structure of the single layer network is the same as the contagion layer of the multiplex network defined above. As for the static network, except that the awareness layer is a static network, its topology structure is totally the same as that of the time-varying network. On each panel, according to different values of the aware probability λ and the forgetting rate δ, we plot four dashed lines for the coupled dynamical processes. The red circle and green diamond line are under the condition when δ = 0.2 and λ = 0.8, while for the red square and green triangle line, δ = 0.8 and λ = 0.2. As for the recovery rate μ, the value is set to be (a) μ = 0.2, (b) μ = 0.4, (c) μ = 0.6, and (d) μ = 0.8. Each line is obtained by averaging 50 independent MC simulations.

Close modal

As for the information spreading, apart from the topology structure of network, the way information spreads is also very important. Since in reality, individuals exhibit herd-like behavior because they are making decisions based on the actions of other individuals, which is called information cascade,39 instead of the classic SIS model, we use a threshold model to simulate the spreading of information. In the threshold model, for an unaware individual, awareness can come from two sources: the ratio between the number of aware neighbors and all its connections which is also known as node degree surpasses the critical value (local threshold α) or it is already infected. In Fig. 5, we illustrate the threshold model defined on the time varying layer.

FIG. 5.

Illustration of the threshold model used to simulate the spreading of information. Only if the ratio between the number of aware neighbors and its degree larger than the threshold value α, the unaware individual can become aware. Here, in the toy model, the threshold is set to be 0.5. The red individuals are aware, while the blue ones are unaware.

FIG. 5.

Illustration of the threshold model used to simulate the spreading of information. Only if the ratio between the number of aware neighbors and its degree larger than the threshold value α, the unaware individual can become aware. Here, in the toy model, the threshold is set to be 0.5. The red individuals are aware, while the blue ones are unaware.

Close modal

Therefore, the multiplex network is made up of two layers with different models on each layer, as described in Ref. 28. For the purpose of calculating the epidemic threshold, we need to change the expression of function ri as follows:

ri(t)=H(αjaji(t)pjA(t)ki(t)),
(9)

where aji(t)=1 if there is a connection between node i and node j at time t, otherwise aji(t)=0. Besides, ki(t) is the degree, the number of its connections, of node i at time t. H(x) is a Heaviside step function, i.e., if x > 0, H(x) = 1, else H(x) = 0. Note that since function qiA(t) and qiU(t) (Eq. (1)) represent the transformation between state S and state I, there is no need to change these two functions. Then after the same derivation, we can also get the numerical result of the epidemic threshold

βc=μΛmax,
(10)

where Λmax is the same as before, i.e., the largest eigenvalue of the matrix H defined following Eq. (6). In Fig. 6, we crosscheck the numerical results with extensive MC simulations on top of the same multiplex network defined above. The results show that the agreement between MMCA method and MC simulations is quite well. Besides, as illustrated in Ref. 28, the two stage effects, which means that there exists a sharp like transition for the epidemic threshold at α0.5, also exists in our model irrespective of the value of δ and μ. Since in our model, the information spreading layer is a time varying network, which is very different from the scale-free network or Erdős-Rényi network, the existence of the two stage effects verifies again that the phenomenon does not rely on the structure of the network, but is a result of the coupled dynamical processes itself. In addition, there is another interesting phenomenon that in Figs. 3 and 6, the agreement between the MMCA approximation and the MC simulations seems to degrade as μ increases and δ decreases. Actually, a large value of μ and a small value of δ can lower the final epidemic size and heighten the epidemic threshold. If we focus on the absolute discrepancy between the MMCA method and MC simulations, it seems that the discrepancy is large on the occasion. However, from a macro point of view, for example, the error rate which equals βMCβMMCAβMC, we find the MMCA method always has a good performance whatever μ and δ is since the error rate locates in the range of 0%–5%.

FIG. 6.

The comparisons of the epidemic threshold between numerical results (MMCA method) and MC simulations as a function of the infected probability β and the local threshold α. The green dotted line corresponds to the epidemic thresholds calculated by the MMCA method, while the heatmap represents the density of the infected nodes ρI obtained by MC simulations under different parameters. From top left to bottom right, the recovery rate μ and the forgetting rate δ are set to be (a) μ = 0.2, δ = 0.8, (b) μ = 0.4, δ = 0.6, (c) μ = 0.6, δ = 0.4, and (d) μ = 0.8, δ = 0.2, respectively. The four phase diagrams are obtained by averaging 100 MC simulations for each point in the grid 50 × 50.

FIG. 6.

The comparisons of the epidemic threshold between numerical results (MMCA method) and MC simulations as a function of the infected probability β and the local threshold α. The green dotted line corresponds to the epidemic thresholds calculated by the MMCA method, while the heatmap represents the density of the infected nodes ρI obtained by MC simulations under different parameters. From top left to bottom right, the recovery rate μ and the forgetting rate δ are set to be (a) μ = 0.2, δ = 0.8, (b) μ = 0.4, δ = 0.6, (c) μ = 0.6, δ = 0.4, and (d) μ = 0.8, δ = 0.2, respectively. The four phase diagrams are obtained by averaging 100 MC simulations for each point in the grid 50 × 50.

Close modal

Moreover, we also compare our model on the multiplex network with other models, including the SIS model on the single layer network and the threshold model on the static multiplex network, as can be seen in Fig. 7. Obviously, as a result of information spreading, the epidemic threshold of our model is larger than that of the single layer case. The final epidemic size is also the smaller one, which means that information spreading can always suppress the spreading of epidemics. The differences of the spreading process between the static multiplex network and the time varying case also crosscheck the findings about time variability of the awareness layer topology, as discussed in Fig. 4.

FIG. 7.

MC simulations for epidemic spreading on multiplex network (circle lines) and single layer network (blue lines), as well as a two-layer static multiplex network (square lines). The topology structure of the single layer network is the same as the contagion layer of the multiplex network defined above. As for the static network, except that the awareness layer is a static network, its topology structure is totally the same as the time-varying network. The recovery rate μ is set as follows: (a) μ=0.2, (b) μ=0.4, (c) μ=0.6, and (d) μ=0.8. In each panel, the dashed lines represent the condition when the forgetting rate δ=0.2 and for the dotted lines, δ=0.8. Besides, the awareness probability is set to be 0.2. Every line is obtained by averaging 100 MC simulations.

FIG. 7.

MC simulations for epidemic spreading on multiplex network (circle lines) and single layer network (blue lines), as well as a two-layer static multiplex network (square lines). The topology structure of the single layer network is the same as the contagion layer of the multiplex network defined above. As for the static network, except that the awareness layer is a static network, its topology structure is totally the same as the time-varying network. The recovery rate μ is set as follows: (a) μ=0.2, (b) μ=0.4, (c) μ=0.6, and (d) μ=0.8. In each panel, the dashed lines represent the condition when the forgetting rate δ=0.2 and for the dotted lines, δ=0.8. Besides, the awareness probability is set to be 0.2. Every line is obtained by averaging 100 MC simulations.

Close modal

From the results above, we can find that the epidemic threshold calculated by the MMCA method has a good agreement with MC simulations. Since according to Eqs. (1)–(3), it is easy for us to obtain the steady density of the infected nodes ρI, here in order to have a better understanding of the MMCA method, we plot ρI as a function of β by means of the MMCA method and MC simulations, as illustrated in Fig. 8. The findings reveal that the agreement between the MMCA method and MC simulations decreases with the increase in μ and the decrease in δ, which is the same as in Figs. 3 and 6.

FIG. 8.

The comparison of the steady density of the infected nodes ρiI between MMCA method (triangles) and MC simulations (circles). Left panel is obtained under the SIS model, while right panel is obtained under the threshold model. The other parameters are set as follows: λ=0.2,α=0.2. Each line of the MC simulations is obtained by averaging 100 independent experiments.

FIG. 8.

The comparison of the steady density of the infected nodes ρiI between MMCA method (triangles) and MC simulations (circles). Left panel is obtained under the SIS model, while right panel is obtained under the threshold model. The other parameters are set as follows: λ=0.2,α=0.2. Each line of the MC simulations is obtained by averaging 100 independent experiments.

Close modal

Since we propose two models to explore the coupled dynamical processes on multiplex network with awareness layer being a time varying network, it is of most interest for us to explore the effects that the varying topology has on different models. As for the activity driven model, m is a critical parameter for the topology structure of the resulting time varying network. Accordingly in the following, we perform many simulations on the two models by setting m to be 4, 7, 10, and 20. Furthermore, in order to study to what extent the coupled dynamical processes is affected by varying topology, a new index named fluctuation ratio is introduced, which is defined as follows:

FR(t)=f(t)f(4)f(4),
(11)

where f(t) is the final epidemic size of the coupled dynamical processes when m is set to be t. It is clear that for both models, with the increase in m, the epidemic threshold and final epidemic size become smaller. The reason is that m represents the number of edges each active node can connect, the larger the m is, the denser is the network. As a result, the spreading process on the dense network is quicker, which makes the nodes become aware easier. This will then suppress the spreading of epidemics. Moreover, the results also show that the threshold model is more susceptible to random changes in the topology of time varying network, since the fluctuation of the final epidemic size, as well as the epidemic threshold, is much larger in this case (Fig. 9).

FIG. 9.

The effects of different time varying topologies on the spreading of epidemics under the SIS model and the threshold model. For each model, we plot four dotted lines through changing the parameter m, namely, 4, 7, 10, and 20, which can produce different time varying networks. The left two panels show the percent of the infected nodes ρI as a function of β. Meanwhile, the right panels correspond to the fluctuation ratio calculated by Eq. (11).

FIG. 9.

The effects of different time varying topologies on the spreading of epidemics under the SIS model and the threshold model. For each model, we plot four dotted lines through changing the parameter m, namely, 4, 7, 10, and 20, which can produce different time varying networks. The left two panels show the percent of the infected nodes ρI as a function of β. Meanwhile, the right panels correspond to the fluctuation ratio calculated by Eq. (11).

Close modal

In addition, as described above, the change of the varying topology of awareness layer can result in different spreading processes for awareness, and then, the epidemic spreading process is also changed. It is also meaningful for us to quantify to what extent the epidemic spreading is affected by the awareness parameter, which is λ for the SIS model or α for the threshold model. Thus, in the following, with the help of MMCA method, we study the evolution of epidemic threshold according to the values of μ and δ. In order to have a better understanding of the effect of λ or α, we plot the phase space of μ and δ, and the color of it corresponds to the value which is the epidemic threshold at λ = 0.8 (α = 0.2) minus the epidemic threshold at λ = 0.2 (α = 0.8), as illustrated in Fig. 10. The results show that the awareness parameter has a more obvious effect on the threshold model since the value of βα=0.2βα=0.8 is much larger than that of the SIS model, especially when δ is small and μ is large. Taking the large fluctuation rate phenomenon studied above into account, we find that the threshold model is “frail” under the configuration of our multiplex network.

FIG. 10.

The phase space of μ and δ with color corresponding to the effect of λ (SIS model, left panel) or α (threshold model, right panel). The settings of the multiplex network are the same as those described in Fig. 3. Importantly, each panel is obtained by two steps: for each pair values of (μ, δ), first, we calculate the epidemic thresholds by letting λ (α) be 0.8 and 0.2, respectively. Then, the final value, which is shown on each panel, equals the epidemic threshold at λ = 0.8 (α = 0.2) minus the epidemic threshold at λ = 0.2 (α = 0.2).

FIG. 10.

The phase space of μ and δ with color corresponding to the effect of λ (SIS model, left panel) or α (threshold model, right panel). The settings of the multiplex network are the same as those described in Fig. 3. Importantly, each panel is obtained by two steps: for each pair values of (μ, δ), first, we calculate the epidemic thresholds by letting λ (α) be 0.8 and 0.2, respectively. Then, the final value, which is shown on each panel, equals the epidemic threshold at λ = 0.8 (α = 0.2) minus the epidemic threshold at λ = 0.2 (α = 0.2).

Close modal

As a summary, in this paper, through considering the interactions among individuals as a time varying network where the awareness diffusion process occurs, we have explored the effects that the spread of awareness has on epidemic spreading. With the help of the multiplex network, we are able to model the interplay between these two kinds of dynamical processes. Afterwards, we build the probability tree for the coupled dynamical processes to reveal the transitions among different states of individuals. Using the MMCA method, the epidemic threshold can be obtained by solving an eigenvalue problem. Regardless of whether the information spreading process is defined as the SIS model or the threshold model, our numerical results of epidemic threshold show high accuracy compared with extensive MC simulations. Specifically, when the information spreading process is the SIS model, the coupled dynamical processes make it difficult for epidemics to outbreak and can also lower the density of the infected individuals. However, if the information spreading process is assumed to be a threshold model, the epidemic threshold exists as a sharp like transition when the local threshold α0.5 though all the epidemic thresholds are smaller than the single layer threshold.

Furthermore, compared with the static multiplex network with the same dynamics, we find that the temporal setting hinders the spreading of awareness which directly affect the epidemic threshold. These results show the importance of taking the information spreading process into account when we try to control the spread of epidemics. From the point of view of information diffusion model, our results show that temporal changes in the topology also hinder the spread of awareness, especially for the threshold model. Since in reality with the burst development of social networks, the way how information spreads is a very complex problem for scientists from many different fields, proposing a proper model to study it is of great significance. Although in this paper, we use the activity driven model to construct the propagation path of information, there are still many other different approaches to realize the goal, for example, bursty model36 and so on. Only if we have a good understanding of the two dynamical processes, which include information spreading and epidemic spreading, it can be possible for us to explore the coupled effects between them. Besides, the results show that changes in the topology of the awareness layer could change the location of the meta-critical point. This is particularly important since a slight delay or advance in the epidemic could have important implications to vaccination strategies.40 Overall, our results give us some useful suggestions about how to model this kind of coupled system and also shed light on new strategies on restraining epidemics.

We thank the anonymous referees for helpful comments and constructive suggestions. This work was partially supported by the NSFC Nos. 11201017 and 11290141, Cultivation Project of NSFC (No. 91130019). Q.G. is also thankful to China Scholarship Council (CSC) (No. 201406020055) for financial support.

1.
R.
Pastor-Satorras
and
A.
Vespignani
,
Phys. Rev. Lett.
86
,
3200
(
2001
).
2.
M. E. J.
Newman
,
Phys. Rev. E
66
,
016128
(
2002
).
3.
Y.
Moreno
,
R.
Pastor-Satorras
, and
A.
Vespignani
,
Eur. Phys. J. B
26
,
521
(
2002
).
4.
Z.
Dezso
and
A.-L.
Barabási
,
Phys. Rev. E
65
,
055103
(
2002
).
5.
D.
Balcan
and
A.
Vespignani
,
Nat. Phys.
7
,
581
(
2011
).
6.
H. E.
Stanley
,
Introduction to Phase Transitions and Critical Phenomena
(
Oxford University Press
,
Oxford
,
1987
).
7.
M.
Barthélemy
,
A.
Barrat
,
R.
Pastor-Satorras
, and
A.
Vespignani
,
Phys. Rev. Lett.
92
,
178701
(
2004
).
8.
J.
Gómez-Gardeñes
,
V.
Latora
,
Y.
Moreno
, and
E.
Profumo
,
Proc. Natl. Acad. Sci. U.S.A.
105
,
1399
(
2008
).
9.
W.
Wang
,
M.
Tang
,
H.
Zhang
, and
Y.
Lai
,
Phys. Rev. E
92
(
1
),
012820
(
2015
).
10.
A.
Motter
and
Y.
Lai
,
Phys. Rev. E
66
(
6
),
065102
(
2002
).
11.
M. E. J.
Newman
,
SIAM Rev.
45
,
167
(
2003
).
12.
S. N.
Dorogovtsev
,
A. V.
Goltsev
, and
J. F. F.
Mendes
,
Rev. Mod. Phys.
80
,
1275
(
2008
).
13.
N. T. J.
Bailey
,
The Mathematical Theory of Infectious Diseases
(
Springer
,
1975
).
14.
A.
Grabowski
and
R. A.
Kosiński
,
Phys. Rev. E
70
,
031908
(
2004
).
15.
M. J.
Keeling
and
K. T. D.
Eames
,
J. R. Soc. Interface
2
,
295
(
2005
).
16.
S.
Funk
,
M.
Salathé
, and
V. A.
Jansen
,
J. R. Soc. Interface
7
,
1247
(
2010
).
17.
M.
Salathé
,
M.
Kazandjieva
,
J. W.
Lee
,
P.
Levis
,
M. W.
Feldman
, and
J. H.
Jones
,
Proc. Natl. Acad. Sci. U.S.A.
107
,
22020
(
2010
).
18.
V. E.
Vergu
,
H.
Busson
, and
P.
Ezanno
,
PloS One
5
,
e9371
(
2010
).
19.
Q.
Wu
,
X.
Fu
,
M.
Small
, and
X.-J.
Xu
,
Chaos
22
,
013101
(
2012
).
20.
S.
Gómez
,
A.
Diaz-Guilera
,
J.
Gómez-Gardeñes
,
C. J.
Pérez-Vicente
,
Y.
Moreno
, and
A.
Arenas
,
Phys. Rev. Lett.
110
,
028701
(
2013
).
21.
P. J.
Mucha
,
T.
Richardson
,
K.
Macon
,
M. A.
Porter
, and
J.-P.
Onnela
,
Science
328
,
876
(
2010
).
22.
M.
Szell
,
R.
Lambiotte
, and
S.
Thurner
,
Proc. Natl. Acad. Sci. U.S.A.
107
,
13636
(
2010
).
23.
Z.
Wang
,
L.
Wang
, and
M.
Perc
,
Phys. Rev. E
89
,
052813
(
2014
).
24.
S.
Boccaletti
,
G.
Bianconi
,
R.
Criado
,
C. I.
del Genio
,
J.
Gómez-Gardeñes
,
M.
Romance
,
I.
Sendiña-Nadal
,
Z.
Wang
, and
M.
Zanin
,
Phys. Rep.
544
,
1
(
2014
).
25.
M.
Kivelä
,
A.
Arenas
,
M.
Barthelemy
,
J. P.
Gleeson
,
Y.
Moreno
, and
M. A.
Porter
,
J. Complex Networks
2
(
3
),
203
(
2014
).
26.
K.-M.
Lee
,
J. Y.
Kim
,
W.-K.
Cho
,
K.-I.
Goh
, and
I.-M.
Kim
,
New J. Phys.
14
,
033027
(
2012
).
27.
C.
Granell
,
S.
Gómez
, and
A.
Arenas
,
Phys. Rev. Lett.
111
,
128701
(
2013
).
28.
Q.
Guo
,
X.
Jiang
,
Y.
Lei
,
M.
Li
,
Y.
Ma
, and
Z.
Zheng
,
Phys. Rev. E
91
,
012822
(
2015
).
29.
C.
Granell
,
S.
Gómez
, and
A.
Arenas
,
Phys. Rev. E
90
,
012808
(
2014
).
30.
S.
Gómez
,
A.
Arenas
,
J.
Borge-Holthoefer
,
S.
Meloni
, and
Y.
Moreno
,
Europhys. Lett.
89
,
38009
(
2010
).
31.
Z.
Wang
,
H.
Zhang
, and
Z.
Wang
,
Chaos
61
,
1
(
2014
).
32.
N.
Perra
,
D.
Balcan
,
B.
Goncalves
, and
A.
Vespignani
,
PloS One
6
,
e23084
(
2011
).
33.
W.
Li
,
S.
Tang
,
W.
Fang
,
Q.
Guo
,
X.
Zhang
, and
Z.
Zheng
,
Phys. Rev. E
92
(
4
),
042810
(
2015
).
34.
N.
Perra
,
B.
Goncalves
,
R.
Pastor-Satorras
, and
A.
Vespignani
,
Sci. Rep.
2
,
469
(
2012
).
35.
M.
Karsai
,
N.
Perra
, and
A.
Vespignani
,
Sci. Rep.
4
,
4001
(
2014
).
36.
E.
Valdano
,
L.
Ferreri
,
C.
Poletto
, and
V.
Colizza
,
Phys. Rev. X
5
,
021005
(
2015
).
37.
Y.
Lei
,
X.
Jiang
,
Q.
Guo
,
Y.
Ma
,
M.
Li
, and
Z.
Zheng
,
Phys. Rev. E
93
,
032308
(
2016
).
38.
P.
Holme
and
J.
Saramäki
,
Phys. Rep.
519
,
97
(
2012
).
39.
J.
Borge-Holthoefer
,
R. A.
Baños
,
S.
González-Bailón
, and
Y.
Moreno
,
J. Complex Networks
1
,
3
(
2013
).
40.
P.
Bajardi
,
C.
Poletto
,
J. J.
Ramasco
,
M.
Tizzoni
,
V.
Colizza
, and
A.
Vespignani
,
PloS One
6
(
1
),
e16591
(
2011
).