The spreading of an infectious disease can trigger human behavior responses to the disease, which in turn plays a crucial role on the spreading of epidemic. In this study, to illustrate the impacts of the human behavioral responses, a new class of individuals, *S ^{F}*, is introduced to the classical susceptible-infected-recovered model. In the model,

*S*state represents that susceptible individuals who take self-initiate protective measures to lower the probability of being infected, and a susceptible individual may go to

^{F}*S*state with a response rate when contacting an infectious neighbor. Via the percolation method, the theoretical formulas for the epidemic threshold as well as the prevalence of epidemic are derived. Our finding indicates that, with the increasing of the response rate, the epidemic threshold is enhanced and the prevalence of epidemic is reduced. The analytical results are also verified by the numerical simulations. In addition, we demonstrate that, because the mean field method neglects the dynamic correlations, a wrong result based on the mean field method is obtained—the epidemic threshold is not related to the response rate, i.e., the additional

^{F}*S*state has no impact on the epidemic threshold.

^{F}Human behavioral responses to the spreading of epidemics have been recognized to have great influence on the epidemic dynamics. Therefore, it is very important to incorporate human behaviors into epidemiological models, which could improve models' utility in reflecting the reality and studying corresponding controlling measures. However, analytically grounded approaches to the problem of interacting effect between epidemic dynamics and human behavioral responses are still lacking so far, and what can help us better understand the impacts of human behavioral responses. In this work, by introducing a new *S*^{F} state into classical susceptible-infected-recovered (*SIR*) model to mimic the situation that, when susceptible individuals are aware of the risk of infection, they may take self-protective measures to lower the probability of being infected. We then derive theoretical formulas based on the percolation method for the epidemic threshold and the prevalence of epidemic. Our results indicate that the introduction of the new *S*^{F} state can significantly enhance the epidemic threshold and reduce the prevalence of epidemic. It is worth mentioning that a wrong result may be obtained when using the traditional mean field method—the epidemic threshold is not altered by the additional *S*^{F} state. The result highlights that if the effects of human behavioral responses are ignored in mathematical modeling, the obtained results cannot really reflect the spreading mechanism of epidemics among human population.

## I. INTRODUCTION

Since network models can well describe the spreading of infectious disease among populations, many efforts have been devoted to studying this field.^{1,2} At first, researchers mainly paid attention to analyze the impact of the network structure on epidemic spreading and the control strategies, for example, how the network topology affects the epidemic threshold and the prevalence of epidemic,^{3–8} or how to design an effective immunization strategy to control the outbreaks of epidemics.^{9,10} In reality, outbreaks of epidemics can trigger spontaneous behavioral responses of individuals to take preventive measures, which in turn alters the epidemic dynamics and affects the disease transmission process.^{11–17} Thus, recently, some studies have made attempts to evaluate the impact of the human behaviors on epidemic dynamics. For instance, Funk *et al.*^{18} studied the impacts of awareness spread on both epidemic threshold and prevalence, and they found that, in a well-mixed population, spread of awareness can reduce the prevalence of epidemic but does not tend to affect the epidemic threshold, yet the epidemic threshold is altered when considering on social networks; Sahneh *et al.* considered a Susceptible-Alter-Infected-Susceptible (*SAIS*) model,^{19} and they found that the way of behavioral response can enhance the epidemic threshold; Meloni *et al.* constructed a meta-population model incorporating several scenarios of self-initiated behavioral changes into the mobility patterns of individuals, and they found that such behavioral changes do not alter the epidemic threshold, but may increase rather than decrease the prevalence of epidemic.^{20} Meanwhile, in Ref. 21, by designing the transmission rate of epidemic as a function of the local infected density or the global infected density, Wu *et al.* investigated the effect of such a behavioral response on the epidemic threshold.

One fact is that the infectious neighbors can infect a susceptible individual, also triggering the awareness of the susceptible individual.^{19,22} In view of this, in Ref. 23, Perra *et al.* introduced a new class of individuals, *S ^{F}*, that represents susceptible people who self-initiate behavioral changes that lead to a reduction in the transmissibility of the infectious disease, into the classical

*SIR*model, and they found that such a model (

*SS*) can induce a rich phase space with multiple epidemic peaks and tipping points. However, the network structure was not incorporated into these models. As we know, when the model is considered within the network based framework, new theoretical tools should be used and new phenomena may be observed. In view of this, we incorporate the network structure into the

^{FIR}*SS*model

^{FIR}^{23}to investigate its spreading dynamics. In the model, when contacting an infectious neighbor, susceptible individuals may be infected (from

*S*state

*I*state) with a transmission rate, or a behavioral response may be triggered (from

*S*state to

*S*state) with a response rate. We provide a theoretical formula for the epidemic threshold as well as the prevalence of epidemic via the percolation method,

^{F}^{5}our results show that the introduction of

*S*class can enhance the epidemic threshold and reduce the prevalence of epidemic. We also demonstrate that a wrong result can be obtained—the introduction of

^{F}*S*class cannot alter the epidemic threshold when using mean field method to such a model. The reasons are presented in Sec. V.

^{F}## II. DESCRIPTIONS OF THE MODEL

For the classical *SIR* model on complex network, where each node on network can be in one of three states: susceptible (*S*), infected (*I*), or recovered (*R*). The transmission rate along each *SI* link is *β*, and an infected node can enter *R* state with a recovery rate *μ*. To reflect the fact that, upon observation of infection, susceptible individuals may adopt protective measures to lower their infection risk, a new class, denoted by *S ^{F}*, is introduced into the original

*SIR*model, we use

*SS*model to denote the modified

^{FIR}*SIR*model in this study. In the model, when an

*S*node contacts an

*I*neighbor, besides the probability of being infected, the

*S*node can go to

*S*state with a response rate $\beta F\u22650$. The transmission rate for the

^{F}*S*nodes is lower; thus, we assume the transmission rate for

^{F}*S*nodes is $\gamma \beta $, where $0\u2264\gamma <1$ is a discount factor.

^{F}The *SS ^{FIR}* model is described by the four following reactions and the associated rates:

Note that the *SS ^{FIR}* model returns to the

*SIR*model once $\beta F=0$, and the

*S*state corresponds to fully vaccinated state when

^{F}*γ*= 0.

## III. THEORETICAL ANALYSIS

In our model, during a sufficiently small time interval $\Delta t$, the transition *rates* of an *SI* edge becoming an *II*, *S ^{FI}*, and

*SR*edge are

*β*,

*β*, and

_{F}*μ*, respectively. As a result, the

*probabilities*of an

*SI*edge becoming an

*II*and

*S*edge are given as $T1=\beta \beta +\beta F+\mu $ and $T2=\beta F\beta +\beta F+\mu $, respectively. Similarly, since the transition rate of an

^{FI}*S*edge becoming an

^{FI}*II*and

*S*during a sufficiently small time interval $\Delta t$ are $\gamma \beta $ and

^{FR}*μ*, the probability of an

*S*edge becoming an

^{FI}*II*edge is $T3=\gamma \beta \gamma \beta +\mu $.

^{24}

To analyze our proposed model, we first define “externally infected neighbor” (EIN) for any node. For node *i*, a neighbor *j* is an EIN, meaning that *j* is infected by its neighbors other than *i* (i.e., *j* is infected even though node *i* is deleted from the networks), which is the basic assumption of the cavity theory in statistical physics. Note that this method is suitable for the networks with negligible number of loops as the network size is sufficiently large^{25}). The probability of neighbor *j* being an EIN of *i* is defined as *θ*; then, for the a node with degree *k*, the probability of having *m* EINs is given as

Let $p(R|m)$ be the probability of *i* being infected when it has number *m* of EIN. To calculate such a probability, we need to recognize that, in our model, an *S* node can become an *I* node through two ways: (a) the *S* node is directly infected; or (b) the *S* node first goes to *S ^{F}* state and then is infected. To facilitate the analysis, we approximately assume that the impacts of

*i*'s infected neighbors on node

*i*happen in a non-overlapping order, i.e., they play roles on node

*i*one by one.

For case (a), the probability of node *i* being infected by the *s*th infected neighbor is given as

Eq. (6) indicates that previous *s* − 1 infected neighbors have not changed the state of node *i* (not become *I* or *S ^{F}* state) before they become

*R*state. Therefore, the probability of

*i*being infected is

For case (b), node *i* should first become *S ^{F}* state, and the probability that the susceptible node

*i*is altered by the

*l*th infected neighbors and becomes

*S*can be written as

^{F}which also indicates that previous *l* − 1 infected neighbors have not changed the state of node *i* before they become *R* state. For the remainder $m\u2212l+1$ infected neighbors (including the infected neighbor who just made *i* go to *S ^{F}* state), they can infect

*i*with probability

As a result, the probability of node *i* first becoming *S ^{F}* state and then going to

*I*state is

The probability $p(R|m)$ is

Then, the EIN probability *θ* is the solution to the self-consistent condition

In Eq. (13), $Q(k)=(k+1)P(k+1)\u3008k\u3009$ is the excess degree distribution, where *P*(*k*) is the degree distribution and $\u3008k\u3009$ is the average degree. The generating function for *Q*(*k*) is given as

There is a trivial solution *θ* = 0 in self-consistent equation (13). To have a non-trivial solution, the following condition must be met:

which implies the epidemic can outbreak when

For the prevalence of epidemic (defined as $R(\u221e)$), we can numerically solve *θ* from self-consistent equation (13), then the formula of $R(\u221e)$ is

In Eq. (17), $G0(.)$ is the generating function of degree distribution *P*(*k*), which is described as

## IV. SIMULATION RESULTS

In this section, we perform an extensive set of Monte Carlo simulations to validate the theoretical predictions in Section III. Here, we carry out simulations on an Erdős-Rényi network (labeled ER network)^{26} with network size *N* = 10 000 and average degree $\u3008k\u3009=10$, and a configuration network generated by an uncorrelated configuration model (UCM).^{27} The configuration network also has *N* = 10 000 nodes and the degree distribution meets $P(k)\u223ck\u22123.0$, whose minimal and maximal degrees are *k _{min}* = 3 and $kmax=N$, respectively.

### A. Results on ER network

Differing from the *SIS* (Susceptible-Infected-Susceptible) model, it is not an easy thing to determine the epidemic threshold for the *SIR* model owing to the non-zero value of *R*. In doing so, in Ref. 28, Shu *et al.* suggested that the variability measure

can well predict the epidemic threshold for the *SIR* model, where *ρ* denotes the prevalence of epidemic in one simulation realization.^{29,30} Δ can be explained as the standard deviation of the epidemic prevalence, and is a standard measure to determine critical point in equilibrium phase on magnetic system.^{31} In our simulations, we have taken at least 1000 independent realizations to predict the epidemic threshold. For convenience, in this study, we set recovery rate $\mu =1.0$.

In Fig. 1, for different response rate *β _{F}*, the value of $R(\u221e)$ (upper panels) and the measure Δ (lower panels) as the functions of the transmission rate

*β*are investigated. As shown in Fig. 1, one can observe that the variability measure can well predict the epidemic threshold for our

*SS*model. As a result, in the following figures, we use this method to determine the epidemic threshold, i.e., the point where the value of Δ is the maximal. Fig. 1 also describes that, no matter $\gamma =0.1$ [see Figs. 1(a) and 1(b)] or $\gamma =0.3$ [see Figs. 1(c) and 1(d)], on the one hand, the epidemic threshold is enhanced as the response rate

^{FIR}*β*is increased. On the other hand, for a fixed value of

_{F}*β*, Figs. 1(a) and 1(c) suggest that the prevalence of epidemic is remarkably reduced when

*β*is increased. The results suggest that, by introducing an additional protective state,

_{F}*S*, to the classical

^{F}*SIR*model, the conclusions are quite different from the previous results which have not incorporated the impacts of human behavioral responses. The result again emphasizes the fact that the spontaneous behavioral responses of individuals to the emergent diseases have vital impacts on the epidemic dynamics. If the behavioral responses are ignored in mathematical modelling, the obtained results cannot really reflect the spreading mechanism of epidemics among human population.

We then compare the theoretical results with the Monte Carlo simulations on ER network in Figs. 2 and 3. Since the degree distribution of an Erdős-Rényi network is $P(k)=e\u2212\u3008k\u3009\u3008k\u3009kk!$, the generating functions meet

According to Ineq. (15), the epidemic threshold *β _{c}* for ER network is determined by the following equation:

Moreover, substituting Eq. (20) into Eqs. (13) and (17), the prevalence of epidemic $R(\u221e)$ can be easily solved.

The comparison of $R(\u221e)$ between the simulations and the theoretical result is plotted in Fig. 2, which indicates that the numerical simulation and the theoretical result are in good agreement. Meanwhile, the epidemic threshold for *β _{c}* obtained from Eq. (21) and from numerical method (i.e., the point where Δ is maximal) is compared in Fig. 3, which also indicates that the epidemic threshold predicated by our method is remarkable agreement with numerical simulations. The result in Fig. 3 also suggests that the epidemic threshold

*β*is increased as the value of

_{c}*γ*is decreased. Importantly, the reduction is more efficient when the response rate

*β*is larger.

_{F}### B. Results on UCM network

Real complex networked systems often possess certain degree of skewness in their degree distributions, typically represented by some scale-free topology. We thus check our model on UCM network with degree distribution $P(k)\u223ck\u22123.0$ to illustrate that our theory can generalize to the networks with heterogenous degree distribution and in the absence of degree-to-degree correlation.

As shown in Figs. 4 and 5, one can see that the analytical results are in good agreement with the numerics. They also indicate that, since increasing *β _{F}* can induce more susceptible individuals go to

*S*state and reducing

^{F}*γ*can lower the risk of

*S*nodes, as a result, both of them can lower the prevalence of epidemic $R(\u221e)$ and increase the epidemic threshold

^{F}*β*.

_{c}## V. MEAN FIELD METHOD FOR THE MODEL

One possible way to describe our proposed $SSFIR$ model is the mean field method, which can be written as

where the factor $\Theta (t)=\u2211k\u2032P(k\u2032|k)Ik\u2032(t)$ represents the probability that any given link points to an infected node. In the absence of any degree correlations $\Theta (t)=1\u3008k\u3009\u2211kkP(k)Ik(t)$.^{32}

Based on the above differential equations, we can obtain that the epidemic threshold $\beta c=\u3008k\u3009\u3008k2\u3009$ (detailed derivation is given in Section VII), which means that the epidemic threshold for our model is not related to the response rate *β _{F}* or the discount factor

*γ*, and which is the same to the epidemic threshold of classical

*SIR*model. The simulation results based on the Eqs. (22)–(25) in Fig. 6 also indicate that, based on mean field method, the epidemic threshold is not altered by different values of

*β*or

_{F}*γ*. However, our previous simulation results based on Monte Carlo method have suggested that the epidemic threshold

*β*is increased when

_{c}*β*is increased or

_{F}*γ*is reduced. That is to say, the conclusion obtained by mean-field method is wrong.

Now let us explain why the mean field method gives a wrong result. It is known that, near the epidemic threshold, the fraction of infected nodes (label *ρ _{I}*) is very small. When using the mean field method, the dynamic correlation is neglected, the probability of

*S*node becoming

*S*is proportional to $O(\rho I)$ since the average fraction of infected nodes among neighborhood equals to $O(\rho I)$. Similarly, the probability of

^{F}*S*node becoming

^{F}*I*node is also proportional to $O(\rho I)$. As a result, the probability of $S\u2192SF\u2192I$ is proportional to $O(\rho I2)$, which leads to the effect of the

*S*is ignored and the epidemic threshold obtained by the mean field method is not related to the value of

^{F}*β*or

_{F}*γ*. In fact, when an

*S*node becomes an

*S*node there must exist at least one infected nodes among the neighborhood of the

^{F}*S*node. More importantly, these infected neighbors may exist for a certain time interval, so the probability of $S\u2192SF\u2192I$ is not proportional to $O(\rho I2)$. However, as deduced in Eq. (9), the dynamics correlation near the epidemic threshold is considered in our above analysis, which can accurately predict the epidemic threshold.

## VI. CONCLUSIONS

To summarize, we have proposed an *SS ^{FIR}* epidemiological model in complex networks, in which the probability of susceptible individuals becoming

*S*state is proportional to the number of infected neighbors, to reflect the fact that individuals are more likely to take protective measures when they find their neighbors are infected. By using theoretical analysis as well as numerical simulations, we found that the prevalence of epidemic is effectively reduced and the epidemic threshold is remarkably increased when the response rate

^{F}*β*is increased or the discount factor

_{F}*γ*is reduced. Moreover, we have demonstrated that the mean field based analysis provides a wrong result: the epidemic threshold is not related to the response rate

*β*or discount factor

_{F}*γ*. The reason is that, near the epidemic threshold, the probability of $S\u2192SF\u2192I$ is a second order infinitesimal since the mean field method ignores the dynamic correlation, which makes the effect of

*S*state to be ignored.

^{F}With the development of technology, information induced awareness or self-protective behaviors can not only diffuse through the contact networks where the diseases spread but also can fast diffuse through many different channels, such as the word of mouth, news media, online social networks, and so on. In view of this, recent well-studied multiplex network theory may be an ideal framework to mimic the interplay of information or related awareness and the epidemic dynamics.^{33–35} Thus, how to generalize our model to multiplex networks and provide theoretical analysis to the model is a challenge in our further work.

## ACKNOWLEDGMENTS

This work was funded by the National Natural Science Foundation of China (Grant Nos. 61473001, 91324002, and 11331009).

### APPENDIX: EPIDEMIC THRESHOLD BASED ON MEAN FIELD METHOD

where $\varphi (t)=\u222b0t\Theta (t\u2032)dt\u2032=1\u3008k\u3009\u2211kkP(k)Rk(t)$.

By using the variation of constants method, $SkF(t)$ is solved as

Then

By letting $\varphi \u221e=limt\u2192\u221e\varphi (t)$, and with the fact that $limt\u2192\u221ed\varphi (t)dt=0$ and $I(t\u2192\u221e)=0$ when the epidemics end, a self-consistent equation can be got from Eq. (A4)

The value $\varphi \u221e=0$ is always a solution. In order to have a non-zero solution, the condition

should be satisfied, which means that the epidemic threshold $\beta c=\u3008k\u3009\u3008k2\u3009$.