Skip to main content

Global behavior of epidemic transmission on heterogeneous networks via two distinct routes

Abstract

In the study of epidemic spreading two natural questions are: whether the spreading of epidemics on heterogenous networks have multiple routes, and whether the spreading of an epidemic is a local or global behavior? In this paper, we answer the above two questions by studying the SIS model on heterogenous networks, and give the global conditions for the endemic state when two distinct routes with uniform rate of infection are considered. The analytical results are also verified by numerical simulations.

Introduction

The dynamical behaviors of epidemic diseases have been studied for quite a long time, SIS and SIR are the two important and fundamental epidemic models [1]. For the SIS epidemic model, each individual can exist in two states: S-susceptible and I-infected. At each time, the susceptible individual which is connected to an infected neighbor will be infected with rate λ. Meanwhile, the infected individuals may be recovered and become susceptible individuals at a rate γ. For the SIR model, there are two different aspects from the SIS model: on one hand, individuals can exist in another state: R-recovered/removed; on the other hand, once an infected individual becomes recovered then the individual cannot be infected again.

In order to explore the mechanism of the evolution of complex networks, in 1999, Barabási and Albert addressed a new model of complex networks: scale-free networks (BA) [2]. In a scale-free network the probability P(k) for any node with k links to other nodes is distributed according to the power law P(k) ~ k-γ. Then researchers found that many real complex systems are scale-free networks, such as the WWW (World Wide Web), the Internet, and so on. These types of networks are considered as heterogenous networks (with a degree distribution exhibiting large fluctuations). Because many epidemic diseases occur in communities exhibiting characteristics consistent with heterogenous networks, many researchers studied the mechanism of the spreading of epidemics on heterogenous networks [3–10] (see also Zhang HF, Fu XC: The spread of epidemics on scale-free networks with nonlinear infectivity. Nonl Anal, Series A, in press, and d'Onofrio A: A note on the global behaviour of the network-based SIS epidemic model. Nonl Anal Real World Appl, in press).

For the mechanism of the spreading of epidemic on complex networks, different researcher gave different explanations [3, 6, 10]. For instance, Pastor-Satorras et al. concluded that the epidemic threshold λ c = 0 for heterogenous networks with sufficiently large size [3], and Zhou et al suggested that the threshold λ c is a constant value, regardless of the size of networks and the degree distribution [6]. Both of the results were obtained just by considering one route of spreading epidemic, and the corresponding results are too special to completely reflect the mechanism of spreading of epidemics. Contrary to the above assumptions, many diseases can be spread in many ways, and they have positive thresholds which are relevant to the degree distribution and the size of networks, e.g, people transmit HIV by having unprotected sex, by receiving infected blood transfusions or, through birth. Furthermore, some epidemics just can prevail in a local place, and some kind of epidemics will globally prevail. In order to better explain the mechanism of the spreading of epidemic on complex networks, and answer whether the spreading of an epidemic is a local or global behavior, we consider two distinct routes of spreading of epidemics on heterogenous networks, and obtain some new results which are not given previously.

The rest of this paper is organized as follows: In Section 2, we study two routes of spreading epidemics on complex networks by a method [11] (and see d'Onofrio A: A note on the global behaviour of the network-based SIS epidemic model. Nonl Anal Real World Appl, in press) which answers explicitly whether the spreading of an epidemic is a local or global behavior. Though the threshold can also be obtained by solving a self-consistency equation [3–5], however, such method did not answer whether the endemic state is locally or globally stable. In order to demonstrate the advantage of the method we used [11] (also d'Onofrio A: A note on the global behaviour of the network-based SIS epidemic model. Nonl Anal Real World Appl, in press), another two routes of spreading epidemics is considered in Section 3, here the threshold cannot be obtained by solving a self-consistency equation. In Section 4, we present numerical simulations to verify our analytical results. And finally, conclusions are given in Section 5.

When the threshold can be obtained by solving a self-consistency equation

As for the mechanism of the spreading of epidemics on heterogeneous networks, R. Pastor-Satorras et al considered that the infective capability of infected individuals is proportional to their degrees, and the factor Θ = ∑ k ′ = 1 k m a x P ( k ′ | k ) I k ′ MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeuiMdeLaeyypa0ZaaabCaeaacqWGqbaucqGGOaakcuWGRbWAgaqbaiabcYha8jabdUgaRjabcMcaPiabdMeajnaaBaaaleaacuWGRbWAgaqbaaqabaaabaGafm4AaSMbauaacqGH9aqpcqaIXaqmaeaacqWGRbWAdaWgaaadbaGaemyBa0MaemyyaeMaemiEaGhabeaaa0GaeyyeIuoaaaa@4374@ (here P (k'|k) denotes the conditional probability that a node with degree k is connected to a node with degree k', k max is the maximal degree of the networks, and I k' is the density of infected individuals with degree k') was used to stand for the probability that an edge emanating from a node of degree k points to an infected nodes. As a result, the epidemic threshold λ c = 0 for heterogenous networks with sufficiently large size [3–5]. However, in [6], Zhou et al argue that Θ may have another form, they let Θ = ∑ k ′ = 1 k m a x U P ( k ′ | k ) I k ′ k ′ MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeuiMdeLaeyypa0ZaaabCaKqbagaadaWcaaqaaiabdwfavjabdcfaqjabcIcaOiqbdUgaRzaafaGaeiiFaWNaem4AaSMaeiykaKIaemysaK0aaSbaaeaacuWGRbWAgaqbaaqabaaabaGafm4AaSMbauaaaaaaleaacuWGRbWAgaqbaiabg2da9iabigdaXaqaaiabdUgaRnaaBaaameaacqWGTbqBcqWGHbqycqWG4baEaeqaaaqdcqGHris5aaaa@46B0@ , just because every infected individual may have the same infective capability U to infect other individuals. And they obtained the epidemic threshold λ c = 1 U MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeq4UdW2aaSbaaSqaaiabdogaJbqabaGccqGH9aqpjuaGdaWcaaqaaiabigdaXaqaaiabdwfavbaaaaa@32D5@ , no matter the individuals' degree and the size of scale-free networks.

Just as we mentioned in above section, both of them have limitations, so we consider that both of these mechanisms are concurrent on complex networks.

For the heterogeneous networks, we also classify individuals into groups according to their degrees, then S k and I k are the density of susceptible individuals and the density of infected individuals with degree k respectively, so S k + I k = 1.

Now we give the dynamical equations for the epidemic on complex networks:

d I k ( t ) d t = − I k + p λ k ( 1 − I k ) Θ 1 + … ( 1 − p β ( 1 − I k ) k Θ 2 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeWabiWaaaqaaKqbaoaalaaabaGaemizaqMaemysaK0aaSbaaeaacqWGRbWAaeqaaiabcIcaOiabdsha0jabcMcaPaqaaiabdsgaKjabdsha0baaaOqaaiabg2da9aqaaiabgkHiTiabdMeajnaaBaaaleaacqWGRbWAaeqaaOGaey4kaSIaemiCaaNaeq4UdWMaem4AaSMaeiikaGIaeGymaeJaeyOeI0IaemysaK0aaSbaaSqaaiabdUgaRbqabaGccqGGPaqkcqqHyoqudaWgaaWcbaGaeGymaedabeaakiabgUcaRiablAcilbqaaaqaaaqaaiabcIcaOiabigdaXiabgkHiTiabdchaWjabek7aIjabcIcaOiabigdaXiabgkHiTiabdMeajnaaBaaaleaacqWGRbWAaeqaaOGaeiykaKIaem4AaSMaeuiMde1aaSbaaSqaaiabikdaYaqabaGccqGGSaalaaaaaa@5C6E@
(1)

for k = 1, 2, ..., k max and where λ, β are the two different infective rates, and the recovery rate is assumed to be unity, the parameter p, (0 ≤ p ≤ 1) gives the different ratio between two routes of spreading epidemic. We suppose the degree distribution is uncorrelated, that is, P ( k ′ | k ) = k ′ P ( k ′ ) 〈 k 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemiuaaLaeiikaGIafm4AaSMbauaacqGG8baFcqWGRbWAcqGGPaqkcqGH9aqpjuaGdaWcaaqaaiqbdUgaRzaafaGaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkaeaacqGHPms4cqWGRbWAcqGHQms8aaaaaa@3F31@ , then we have

Θ 1 = ∑ k ′ = 1 k m a x k ′ P ( k ′ ) I k ′ 〈 k 〉 Θ 2 = ∑ k ′ = 1 k m a x U P ( k ′ ) I k ′ 〈 k 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeqabiqaaaqaaiabfI5arnaaBaaaleaacqaIXaqmaeqaaOGaeyypa0ZaaabCaKqbagaadaWcaaqaaiqbdUgaRzaafaGaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkcqWGjbqsdaWgaaqaaiqbdUgaRzaafaaabeaaaeaacqGHPms4cqWGRbWAcqGHQms8aaaaleaacuWGRbWAgaqbaiabg2da9iabigdaXaqaaiabdUgaRnaaBaaameaacqWGTbqBcqWGHbqycqWG4baEaeqaaaqdcqGHris5aaGcbaGaeuiMde1aaSbaaSqaaiabikdaYaqabaGccqGH9aqpdaaeWbqcfayaamaalaaabaGaemyvauLaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkcqWGjbqsdaWgaaqaaiqbdUgaRzaafaaabeaaaeaacqGHPms4cqWGRbWAcqGHQms8aaaaleaacuWGRbWAgaqbaiabg2da9iabigdaXaqaaiabdUgaRnaaBaaameaacqWGTbqBcqWGHbqycqWG4baEaeqaaaqdcqGHris5aaaaaaa@65A6@
(2)

where 〈 k 〉 = ∑ k ′ = 1 k m a x k ′ P ( k ′ ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeyykJeUaem4AaSMaeyOkJeVaeyypa0ZaaabCaeaacuWGRbWAgaqbaiabdcfaqjabcIcaOiqbdUgaRzaafaGaeiykaKcaleaacuWGRbWAgaqbaiabg2da9iabigdaXaqaaiabdUgaRnaaBaaameaacqWGTbqBcqWGHbqycqWG4baEaeqaaaqdcqGHris5aaaa@42C4@ . Then Eq(1) can be written in a compact form:

d I k ( t ) d t = − I k + … ( 1 − I k ) k ∑ k ′ = 1 k m a x ( p λ k ′ + ( 1 − p ) U β ) P ( k ′ ) I k ′ 〈 k 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeaabiWaaaqaaKqbaoaalaaabaGaemizaqMaemysaK0aaSbaaeaacqWGRbWAaeqaaiabcIcaOiabdsha0jabcMcaPaqaaiabdsgaKjabdsha0baaaOqaaiabg2da9aqaaiabgkHiTiabdMeajnaaBaaaleaacqWGRbWAaeqaaOGaey4kaSIaeSOjGSeabaaabaaabaGaeiikaGIaeGymaeJaeyOeI0IaemysaK0aaSbaaSqaaiabdUgaRbqabaGccqGGPaqkcqWGRbWAdaaeWbqcfayaamaalaaabaGaeiikaGIaemiCaaNaeq4UdWMafm4AaSMbauaacqGHRaWkcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGvbqvcqaHYoGycqGGPaqkcqWGqbaucqGGOaakcuWGRbWAgaqbaiabcMcaPiabdMeajnaaBaaabaGafm4AaSMbauaaaeqaaaqaaiabgMYiHlabdUgaRjabgQYiXdaaaSqaaiqbdUgaRzaafaGaeyypa0JaeGymaedabaGaem4AaS2aaSbaaWqaaiabd2gaTjabdggaHjabdIha4bqabaaaniabggHiLdaaaaaa@6B9E@
(3)

(where k = 1, 2, ..., k max ). By letting I = [I1, I2, ..., I = [ I 1 , I 2 , ⋯ , I k m a x ] T MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemysaKKaeyypa0Jaei4waSLaemysaK0aaSbaaSqaaiabigdaXaqabaGccqGGSaalcqWGjbqsdaWgaaWcbaGaeGOmaidabeaakiabcYcaSiabl+UimjabcYcaSiabdMeajnaaBaaaleaacqWGRbWAdaWgaaadbaacbiGae8xBa0Mae8xyaeMae8hEaGhabeaaaSqabaGccqGGDbqxdaahaaWcbeqaaiabdsfaubaaaaa@41F4@ , Eq(3) can be rewritten as a vector form:

d I ( t ) d t = A I + N ( I , t ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqcfa4aaSaaaeaacqWGKbazcqWGjbqscqGGOaakcqWG0baDcqGGPaqkaeaacqWGKbazcqWG0baDaaGaeyypa0JaemyqaeKaemysaKKaey4kaSIaemOta4KaeiikaGIaemysaKKaeiilaWIaemiDaqNaeiykaKcaaa@3F63@
(4)

where AI is the linear part of I, and N(I, t) is the nonlinear part of I, and

A k k ′ = − δ k k ′ + k ( p λ k ′ + U ( 1 − p ) β ) P ( k ′ ) 〈 k 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemyqae0aaSbaaSqaaiabdUgaRjqbdUgaRzaafaaabeaakiabg2da9iabgkHiTiabes7aKnaaBaaaleaacqWGRbWAcuWGRbWAgaqbaaqabaGccqGHRaWkjuaGdaWcaaqaaiabdUgaRjabcIcaOiabdchaWjabeU7aSjqbdUgaRzaafaGaey4kaSIaemyvauLaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeiykaKIaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkaeaacqGHPms4cqWGRbWAcqGHQms8aaaaaa@51B5@
(5)

(k, k' = 1, 2, ..., k max ) where δ kk' = 1 if k = k', or δ kk' = 0 otherwise.

N k = − I k k ∑ k ′ = 1 k m a x ( p λ k ′ + U ( 1 − p ) β ) P ( k ′ ) I k ′ 〈 . k 〉 < 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeaabiWaaaqaaiabd6eaonaaBaaaleaacqWGRbWAaeqaaaGcbaGaeyypa0dabaGaeyOeI0IaemysaK0aaSbaaSqaaiabdUgaRbqabaGccqWGRbWAdaaeWbqcfayaamaalaaabaGaeiikaGIaemiCaaNaeq4UdWMafm4AaSMbauaacqGHRaWkcqWGvbqvcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqaHYoGycqGGPaqkcqWGqbaucqGGOaakcuWGRbWAgaqbaiabcMcaPiabdMeajnaaBaaabaGafm4AaSMbauaaaeqaaaqaaiabgMYiHlabc6caUiabdUgaRjabgQYiXdaaaSqaaiqbdUgaRzaafaGaeyypa0JaeGymaedabaGaem4AaS2aaSbaaWqaaiabd2gaTjabdggaHjabdIha4bqabaaaniabggHiLdaakeaaaeaacqGH8aapaeaacqaIWaamaaaaaa@5E7C@
(6)

(k = 1, 2, ..., k max ). By some computation, we can find that, for the matrix A, there are (k max - 1) eigenvalues equal to -1.

In order to find the last eigenvalue of A, we let

V = [λ, 2λ, ..., Nλ]T

and write matrix A as

A = − E k m a x + … 1 〈 k 〉 [ ( p λ + U ( 1 − p ) β ) P ( 1 ) V ( 2 p λ + U ( 1 − p ) β ) P ( 2 ) V ⋮ ( N p λ + U ( 1 − p ) β ) P ( N ) V ] T MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeaabiWaaaqaaiabdgeabbqaaiabg2da9aqaaiabgkHiTiabdweafnaaBaaaleaacqWGRbWAdaWgaaadbaGaemyBa0MaemyyaeMaemiEaGhabeaaaSqabaGccqGHRaWkcqWIMaYsaeaaaeaaaeaajuaGdaWcaaqaaiabigdaXaqaaiabgMYiHlabdUgaRjabgQYiXdaakmaadmaabaqbaeqabqqaaaaabaGaeiikaGIaemiCaaNaeq4UdWMaey4kaSIaemyvauLaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeiykaKIaemiuaaLaeiikaGIaeGymaeJaeiykaKIaemOvayfabaGaeiikaGIaeGOmaiJaemiCaaNaeq4UdWMaey4kaSIaemyvauLaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeiykaKIaemiuaaLaeiikaGIaeGOmaiJaeiykaKIaemOvayfabaGaeSO7I0eabaGaeiikaGIaemOta4KaemiCaaNaeq4UdWMaey4kaSIaemyvauLaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeiykaKIaemiuaaLaeiikaGIaemOta4KaeiykaKIaemOvayfaaaGaay5waiaaw2faamaaCaaaleqabaGaemivaqfaaaaaaaa@7DCC@
(8)

where E k m a x MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemyrau0aaSbaaSqaaiabdUgaRnaaBaaameaacqWGTbqBcqWGHbqycqWG4baEaeqaaaWcbeaaaaa@32D2@ is an identity matrix, then we have

A V = ( − 1 + … 1 〈 k 〉 ∑ k ′ = 1 k m a x ( k ′ 2 p λ + U k ′ ( 1 − p ) β ) P ( k ′ ) ) V = ( − 1 + p λ 〈 k 2 〉 〈 k 〉 + ( 1 − p ) U β ) V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeaabmWaaaqaaiabdgeabjabdAfawbqaaiabg2da9aqaaiabcIcaOiabgkHiTiabigdaXiabgUcaRiablAcilbqaaaqaaaqaaKqbaoaalaaabaGaeGymaedabaGaeyykJeUaem4AaSMaeyOkJepaaOWaaabCaeaacqGGOaakcuWGRbWAgaqbamaaCaaaleqabaGaeGOmaidaaOGaemiCaaNaeq4UdWMaey4kaSIaemyvauLafm4AaSMbauaacqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqaHYoGycqGGPaqkcqWGqbaucqGGOaakcuWGRbWAgaqbaiabcMcaPiabcMcaPiabdAfawbWcbaGafm4AaSMbauaacqGH9aqpcqaIXaqmaeaacqWGRbWAdaWgaaadbaGaemyBa0MaemyyaeMaemiEaGhabeaaa0GaeyyeIuoaaOqaaaqaaiabg2da9aqaaiabcIcaOiabgkHiTiabigdaXiabgUcaRiabdchaWjabeU7aSLqbaoaalaaabaGaeyykJeUaem4AaS2aaWbaaeqabaGaeGOmaidaaiabgQYiXdqaaiabgMYiHlabdUgaRjabgQYiXdaakiabgUcaRiabcIcaOiabigdaXiabgkHiTiabdchaWjabcMcaPiabdwfavjabek7aIjabcMcaPiabdAfawbaaaaa@7BDA@
(9)

From Eq (9), it follows that the k max 's eigenvalue of the matrix A is

u = − 1 + p λ 〈 k 2 〉 〈 k 〉 + ( 1 − p ) U β MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemyDauNaeyypa0JaeyOeI0IaeGymaeJaey4kaSIaemiCaaNaeq4UdWwcfa4aaSaaaeaacqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaey4kaSIaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaemyvauLaeqOSdigaaa@48A6@
(10)

If the solution I = 0 of the Eq(3) is stable, all eigenvalues of the matrix A must be non-positive, that is, u ≤ 0.

So we have the following conclusion:

If

p λ 〈 k 2 〉 〈 k 〉 + ( 1 − p ) U β ≤ 1 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemiCaaNaeq4UdWwcfa4aaSaaaeaacqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaey4kaSIaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaemyvauLaeqOSdiMaeyizImQaeGymaedaaa@4613@
(11)

then the solution I = 0 of the Eq(3) is globally asymptotically stable, otherwise the unique endemic solution I = [ I 1 , I 2 , ... , I k m a x ] > 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemysaKKaeyypa0Jaei4waSLaemysaK0aaSbaaSqaaiabigdaXaqabaGccqGGSaalcqWGjbqsdaWgaaWcbaGaeGOmaidabeaakiabcYcaSiabc6caUiabc6caUiabc6caUiabcYcaSiabdMeajnaaBaaaleaacqWGRbWAdaWgaaadbaGaemyBa0MaemyyaeMaemiEaGhabeaaaSqabaGccqGGDbqxcqGH+aGpcqaIWaamaaa@434B@ ( 0 is a zero vector) is globally asymptotically stable, except when there are no infected individuals in the networks at the initial time.

From the inequality (11), we can find the thresholds for the outbreak of epidemics on complex networks with the two routes of spreading epidemic mentioned above, and at the same time we can answer whether the endemic state is globally stable.

When p = 1, that is, there is only one way of spreading of epidemic on complex networks, as discussed in [3], then we obtain the threshold for λ > 〈 k 〉 〈 k 2 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeq4UdWMaeyOpa4tcfa4aaSaaaeaacqGHPms4cqWGRbWAcqGHQms8aeaacqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepaaaaa@3A07@ , the same as the result given in [3–5].

When p = 0, then we obtain the threshold for β > 1 U MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeqOSdiMaeyOpa4tcfa4aaSaaaeaacqaIXaqmaeaacqWGvbqvaaaaaa@313F@ , the same as the result given in [6].

When 0 <p < 1, i.e., there are two routes of spreading of epidemic on complex networks, our result suggest that the threshold for the outbreak of epidemic is positive, which is relevant to the ratio of two routes of spreading of epidemic, the degree distribution, and the size of the networks. So the threshold we obtain is neither zero nor just a constant given by previous authors.

When threshold cannot be obtained by solving a self-consistency equation

In the above section, the threshold for the case we considered can also be obtained by solving a self-consistency equation, although this does not explicitly determine whether the equilibrium is globally stable. In this section, in order to demonstrate the advantage of the method we used in the above section, we consider another case where the threshold cannot be obtained by solving a self-consistency equation. Moreover, such case may exist in the real networks.

We suppose that there are also two routes of the spreading epidemics on complex networks, one route is given by Θ1, another route of the spreading of epidemic is the standard SIS model [1].

The dynamical equations are

d I k ( t ) d t = − I k + … λ p ( 1 − I k ) k 〈 k 〉 ∑ k ′ = 1 k m a x k ′ P ( k ′ ) I k ′ + … ( 1 − I k ) ( 1 − p ) β ( k ) ∑ k ′ = 1 k m a x P ( k ′ ) I k ′ MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeaabmWaaaqcfayaamaalaaabaGaemizaqMaemysaK0aaSbaaeaacqWGRbWAaeqaaiabcIcaOiabdsha0jabcMcaPaqaaiabdsgaKjabdsha0baaaOqaaiabg2da9aqaaiabgkHiTiabdMeajnaaBaaaleaacqWGRbWAaeqaaOGaey4kaSIaeSOjGSeabaaabaaabaGaeq4UdWMaemiCaaNaeiikaGIaeGymaeJaeyOeI0IaemysaK0aaSbaaSqaaiabdUgaRbqabaGccqGGPaqkjuaGdaWcaaqaaiabdUgaRbqaaiabgMYiHlabdUgaRjabgQYiXdaakmaaqahabaGafm4AaSMbauaacqWGqbaucqGGOaakcuWGRbWAgaqbaiabcMcaPiabdMeajnaaBaaaleaacuWGRbWAgaqbaaqabaGccqGHRaWkcqWIMaYsaSqaaiqbdUgaRzaafaGaeyypa0JaeGymaedabaGaem4AaS2aaSbaaWqaaiabd2gaTjabdggaHjabdIha4bqabaaaniabggHiLdaakeaaaeaaaeaacqGGOaakcqaIXaqmcqGHsislcqWGjbqsdaWgaaWcbaGaem4AaSgabeaakiabcMcaPiabcIcaOiabigdaXiabgkHiTiabdchaWjabcMcaPiabek7aIjabcIcaOiabdUgaRjabcMcaPmaaqahabaGaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkcqWGjbqsdaWgaaWcbaGafm4AaSMbauaaaeqaaaqaaiqbdUgaRzaafaGaeyypa0JaeGymaedabaGaem4AaS2aaSbaaWqaaiabd2gaTjabdggaHjabdIha4bqabaaaniabggHiLdaaaaaa@859B@
(12)

(k = 1, ..., k max ). Here, we should distinguish the second term of the right-hand side of Eq(12) from the second term of the right-hand side of Eq(3). In Eq(3), Zhou et al considered that the susceptible individuals which may be infected are proportional to their degrees k though they are irrelevant to the infected individuals, however, in Eq(12) we suppose that all of the susceptible individuals which may be infected are proportional to a constant rate, i.e. (k).

By using the same method as in Section 2, we can rewrite Eq(12) in the following vector form,

d I ( t ) d t = A I + N ( I , t ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqcfa4aaSaaaeaacqWGKbazcqWGjbqscqGGOaakcqWG0baDcqGGPaqkaeaacqWGKbazcqWG0baDaaGccqGH9aqpcqWGbbqqcqWGjbqscqGHRaWkcqWGobGtcqGGOaakcqWGjbqscqGGSaalcqWG0baDcqGGPaqkaaa@3F6D@
(13)

and

A k k ′ = − δ k k ′ + k p λ k ′ P ( k ′ ) 〈 k 〉 + 〈 k 〉 ( 1 − p ) β P ( k ′ ) MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemyqae0aaSbaaSqaaiabdUgaRjqbdUgaRzaafaaabeaakiabg2da9iabgkHiTiabes7aKnaaBaaaleaacqWGRbWAcuWGRbWAgaqbaaqabaGccqGHRaWkjuaGdaWcaaqaaiabdUgaRjabdchaWjabeU7aSjqbdUgaRzaafaGaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkaeaacqGHPms4cqWGRbWAcqGHQms8aaGccqGHRaWkcqGHPms4cqWGRbWAcqGHQms8cqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqaHYoGycqWGqbaucqGGOaakcuWGRbWAgaqbaiabcMcaPaaa@5802@
(14)
N k = − I k ∑ k ′ = 1 k m a x [ k p λ k ′ 〈 k 〉 + 〈 k 〉 ( 1 − p ) β ] P ( k ′ ) I k ′ < 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeaabiWaaaqaaiabd6eaonaaBaaaleaacqWGRbWAaeqaaaGcbaGaeyypa0dabaGaeyOeI0IaemysaK0aaSbaaSqaaiabdUgaRbqabaGcdaaeWbqcfayaaiabcUfaBnaalaaabaGaem4AaSMaemiCaaNaeq4UdWMafm4AaSMbauaaaeaacqGHPms4cqWGRbWAcqGHQms8aaGaey4kaSIaeyykJeUaem4AaSMaeyOkJeVaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeiyxa0LaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkcqWGjbqsdaWgaaqaaiqbdUgaRzaafaaabeaaaSqaaiqbdUgaRzaafaGaeyypa0JaeGymaedabaGaem4AaS2aaSbaaWqaaiabd2gaTjabdggaHjabdIha4bqabaaaniabggHiLdaakeaaaeaacqGH8aapaeaacqaIWaamaaaaaa@6215@
(15)

(k, k' = 1, 2, ..., k max ).

By solving the matrix A, we can find that there are k max - 2 eigenvalues equal to -1, and the other two eigenvalues u are given by the following equation

u 2 − [ λ p + ( 1 − p ) β 〈 k 〉 + λ p 〈 k 〉 T 1 − 2 ] u + … { ( 1 − p ) p β λ T 1 − λ p − ( 1 − p ) β 〈 k 〉 − … λ p 〈 k 〉 T 1 + 1 − ( 1 − p ) p λ β 〈 k 〉 T 2 } = 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeqabmqaaaqaaiabdwha1naaCaaaleqabaGaeGOmaidaaOGaeyOeI0Iaei4waSLaeq4UdWMaemiCaaNaey4kaSIaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeyykJeUaem4AaSMaeyOkJeVaey4kaSscfa4aaSaaaeaacqaH7oaBcqWGWbaCaeaacqGHPms4cqWGRbWAcqGHQms8aaGaemivaq1aaSbaaeaacqaIXaqmaeqaaOGaeyOeI0IaeGOmaiJaeiyxa0LaemyDauNaey4kaSIaeSOjGSeabaGaei4EaSNaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaemiCaaNaeqOSdiMaeq4UdWMaemivaq1aaSbaaSqaaiabigdaXaqabaGccqGHsislcqaH7oaBcqWGWbaCcqGHsislcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqaHYoGycqGHPms4cqWGRbWAcqGHQms8cqGHsislcqWIMaYsaeGabaqhbiaaxMaajuaGdaWcaaqaaiabeU7aSjabdchaWbqaaiabgMYiHlabdUgaRjabgQYiXdaacqWGubavdaWgaaqaaiabigdaXaqabaGccqGHRaWkcqaIXaqmcqGHsislcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGWbaCcqaH7oaBcqaHYoGycqGHPms4cqWGRbWAcqGHQms8cqWGubavdaWgaaWcbaGaeGOmaidabeaakiabc2ha9jabg2da9iabicdaWaaaaaa@95B4@
(16)

where

T 1 = ∑ k ′ = 2 k m a x k ′ ( 1 − k ′ ) P ( k ′ ) = 〈 k 2 〉 − 〈 k 〉 . T 2 = ∑ k ′ = 2 k m a x ( 1 − k ′ ) P ( k ′ ) = 〈 k 〉 − 1 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeqabqWaaaaabaGaemivaq1aaSbaaSqaaiabigdaXaqabaaakeaacqGH9aqpaeaadaaeWbqaaiqbdUgaRzaafaGaeiikaGIaeGymaeJaeyOeI0Iafm4AaSMbauaacqGGPaqkcqWGqbaucqGGOaakcuWGRbWAgaqbaiabcMcaPaWcbaGafm4AaSMbauaacqGH9aqpcqaIYaGmaeaacqWGRbWAdaWgaaadbaGaemyBa0MaemyyaeMaemiEaGhabeaaa0GaeyyeIuoaaOqaaaqaaiabg2da9aqaaiabgMYiHlabdUgaRnaaCaaaleqabaGaeGOmaidaaOGaeyOkJeVaeyOeI0IaeyykJeUaem4AaSMaeyOkJeVaeiOla4cabaGaemivaq1aaSbaaSqaaiabikdaYaqabaaakeaacqGH9aqpaeaadaaeWbqaaiabcIcaOiabigdaXiabgkHiTiqbdUgaRzaafaGaeiykaKIaemiuaaLaeiikaGIafm4AaSMbauaacqGGPaqkaSqaaiqbdUgaRzaafaGaeyypa0JaeGOmaidabaGaem4AaS2aaSbaaWqaaiabd2gaTjabdggaHjabdIha4bqabaaaniabggHiLdaakeaaaeaacqGH9aqpaeaacqGHPms4cqWGRbWAcqGHQms8cqGHsislcqaIXaqmaaaaaa@7332@
(17)

From Eq(16), (17), we have

u 2 − [ ( 1 − p ) β 〈 k 〉 + λ p 〈 k 2 〉 〈 k 〉 − 2 ] u + … { ( 1 − p ) p β λ 〈 k 2 〉 − ( 1 − p ) β 〈 k 〉 − … λ p 〈 k 2 〉 〈 k 〉 + 1 − ( 1 − p ) p λ β 〈 k 〉 2 } = 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeqabmqaaaqaaiabdwha1naaCaaaleqabaGaeGOmaidaaOGaeyOeI0Iaei4waSLaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeyykJeUaem4AaSMaeyOkJeVaey4kaSscfa4aaSaaaeaacqaH7oaBcqWGWbaCcqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaeyOeI0IaeGOmaiJaeiyxa0LaemyDauNaey4kaSIaeSOjGSeabaGaei4EaSNaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaemiCaaNaeqOSdiMaeq4UdWMaeyykJeUaem4AaS2aaWbaaSqabeaacqaIYaGmaaGccqGHQms8cqGHsislcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqaHYoGycqGHPms4cqWGRbWAcqGHQms8cqGHsislcqWIMaYsaeGabaqhbiaaxMaajuaGdaWcaaqaaiabeU7aSjabdchaWjabgMYiHlabdUgaRnaaCaaabeqaaiabikdaYaaacqGHQms8aeaacqGHPms4cqWGRbWAcqGHQms8aaGccqGHRaWkcqaIXaqmcqGHsislcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGWbaCcqaH7oaBcqaHYoGycqGHPms4cqWGRbWAcqGHQms8daahaaWcbeqaaiabikdaYaaakiabc2ha9jabg2da9iabicdaWaaaaaa@9797@
(18)

In order to obtain both of the negative eigenvalues, by solving Eq(18) we have

If

1 ≥ 1 2 [ λ p 〈 k 2 〉 〈 k 〉 + ( 1 − p ) β 〈 k 〉 ] + … 1 2 { [ λ p 〈 k 2 〉 〈 k 〉 − ( 1 − p ) β 〈 k 〉 ] 2 + 4 ( 1 − p ) p λ β 〈 k 〉 2 } 1 / 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeWabiqaaaqaaiabigdaXiabgwMiZMqbaoaalaaabaGaeGymaedabaGaeGOmaidaaOGaei4waSvcfa4aaSaaaeaacqaH7oaBcqWGWbaCcqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaey4kaSIaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeyykJeUaem4AaSMaeyOkJeVaeiyxa0Laey4kaSIaeSOjGSeabaqcfaOaaGPaVlaaykW7daWcaaqaaiabigdaXaqaaiabikdaYaaakmaacmaabaGaei4waSvcfa4aaSaaaeaacqaH7oaBcqWGWbaCcqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaeyOeI0IaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeqOSdiMaeyykJeUaem4AaSMaeyOkJeVaeiyxa01aaWbaaSqabeaacqaIYaGmaaGccqGHRaWkcqaI0aancqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGWbaCcqaH7oaBcqaHYoGycqGHPms4cqWGRbWAcqGHQms8daahaaWcbeqaaiabikdaYaaaaOGaay5Eaiaaw2haamaaCaaaleqabaGaeGymaeJaei4la8IaeGOmaidaaaaaaaa@8BFD@
(19)

then the solution I = 0 of the Eq(15) is globally asymptotically stable, otherwise the unique endemic solution I = [ I 1 , I 2 , ⋯ , I k m a x ] ≠ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemysaKKaeyypa0Jaei4waSLaemysaK0aaSbaaSqaaiabigdaXaqabaGccqGGSaalcqWGjbqsdaWgaaWcbaGaeGOmaidabeaakiabcYcaSiabl+UimjabcYcaSiabdMeajnaaBaaaleaacqWGRbWAdaWgaaadbaGaemyBa0MaemyyaeMaemiEaGhabeaaaSqabaGccqGGDbqxcqGHGjsUcqaIWaamaaa@434C@ ( 0 is a zero vector) is globally asymptotically stable. Except when there are no infected individuals in the networks at the initial time.

When p = 1, then we obtain the threshold for λ > 〈 k 〉 〈 k 2 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeq4UdWMaeyOpa4tcfa4aaSaaaeaacqGHPms4cqWGRbWAcqGHQms8aeaacqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepaaaaa@3A07@ , the same as the result given in [3–5]. When p = 0, then we obtain the threshold for β > 1 〈 k 〉 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdbba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeqOSdiMaeyOpa4tcfa4aaSaaaeaacqaIXaqmaeaacqGHPms4cqWGRbWAcqGHQms8aaaaaa@34EE@ , the same as the result on homogenous networks.

If 0 <p < 1, we can find the threshold for this two routes of spreading epidemic is not just the linear combination of threshold λ and β as shown in inequality (11), which manifest that some complex behaviors may come forth.

We remark here that we can also consider more routes' spreading of epidemics on complex networks, by using the similar method, although here in this paper we have only considered two routes of spreading epidemics.

Numerical simulations

In this section, we present numerical simulations to support the results obtained in previous sections. Our simulations are based on the BA network with P(k) = k-γ, γ = 3, N = 2000, and ⟨k⟩ = 6, ⟨k2⟩ = 98.

In Figure 1, our simulations verify the inequality(11). In order to simultaneously consider two infective rates λ, β, in Figure 1(a), we study the ratio of λ/β with p changing, then inequality (11) can be written as

Figure 1
figure 1

Verification of inequalities (20) and (21). Figure 1 (a) Simulations show that λ/β is in accordance with inequality (20) by changing p. Figure 1 (b) Simulations show that β/λ is in accordance with inequality (21) by changing p.

p ( λ β ) 〈 k 2 〉 〈 k 〉 + ( 1 − p ) U ≤ 1 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemiCaaNaeiikaGscfa4aaSaaaeaacqaH7oaBaeaacqaHYoGyaaGccqGGPaqkjuaGdaWcaaqaaiabgMYiHlabdUgaRnaaCaaabeqaaiabikdaYaaacqGHQms8aeaacqGHPms4cqWGRbWAcqGHQms8aaGccqGHRaWkcqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGvbqvcqGHKjYOcqaIXaqmaaa@486D@
(20)

In Figure 1(b), we study the ratio of β/λ, with p changing, then inequality (11) can be written as

p 〈 k 2 〉 〈 k 〉 + ( 1 − p ) U ( β λ ) ≤ 1 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaemiCaaxcfa4aaSaaaeaacqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaey4kaSIaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaemyvauLaeiikaGscfa4aaSaaaeaacqaHYoGyaeaacqaH7oaBaaGccqGGPaqkcqGHKjYOcqaIXaqmaaa@486D@
(21)

In Figure 2, our simulations verify the inequality (19). Just as in Figure 1, we simultaneously consider two infective rates λ, β, in Figure 2(a), we study the ratio of λ/β with p changing, then inequality (19) can be written as

Figure 2
figure 2

Verification of inequalities (22) and (23). Figure 2 (a) Simulations show that λ/β is in accordance with inequality (22) by changing p. Figure 2 (b) Simulations show that β/λ is in accordance with inequality (23) by changing p.

1 ≥ 1 2 [ λ β p 〈 k 2 〉 〈 k 〉 + ( 1 − p ) 〈 k 〉 ] + … 1 2 { [ λ β p 〈 k 2 〉 〈 k 〉 − ( 1 − p ) 〈 k 〉 ] 2 + λ β 4 ( 1 − p ) p 〈 k 〉 2 } 1 / 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeWabiqaaaqaaiabigdaXiabgwMiZMqbaoaalaaabaGaeGymaedabaGaeGOmaidaaOGaei4waSvcfa4aaSaaaeaacqaH7oaBaeaacqaHYoGyaaWaaSaaaeaacqWGWbaCcqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaey4kaSIaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeyykJeUaem4AaSMaeyOkJeVaeiyxa0Laey4kaSIaeSOjGSeabaqcfaOaaGPaVlaaykW7daWcaaqaaiabigdaXaqaaiabikdaYaaakmaacmaabaGaei4waSvcfa4aaSaaaeaacqaH7oaBaeaacqaHYoGyaaWaaSaaaeaacqWGWbaCcqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaeyOeI0IaeiikaGIaeGymaeJaeyOeI0IaemiCaaNaeiykaKIaeyykJeUaem4AaSMaeyOkJeVaeiyxa01aaWbaaSqabeaacqaIYaGmaaGccqGHRaWkdaWcaaqaaiabeU7aSbqaaiabek7aIbaacqaI0aancqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGWbaCcqGHPms4cqWGRbWAcqGHQms8daahaaWcbeqaaiabikdaYaaaaOGaay5Eaiaaw2haamaaCaaaleqabaGaeGymaeJaei4la8IaeGOmaidaaaaaaaa@8C2D@
(22)

In Figure 2(b), we study the ratio of β/λ, with p changing, then inequality (19) can be written as

1 ≥ 1 2 [ p 〈 k 2 〉 〈 k 〉 + β λ ( 1 − p ) 〈 k 〉 ] + … 1 2 { [ p 〈 k 2 〉 〈 k 〉 − β λ ( 1 − p ) 〈 k 〉 ] 2 + β λ 4 ( 1 − p ) p 〈 k 〉 2 } 1 / 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qqaqFr0xc9vqFj0dXdbba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaqbaeWabiqaaaqaaiabigdaXiabgwMiZMqbaoaalaaabaGaeGymaedabaGaeGOmaidaaOGaei4waSvcfa4aaSaaaeaacqWGWbaCcqGHPms4cqWGRbWAdaahaaqabeaacqaIYaGmaaGaeyOkJepabaGaeyykJeUaem4AaSMaeyOkJepaaOGaey4kaSscfa4aaSaaaeaacqaHYoGyaeaacqaH7oaBaaGccqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqGHPms4cqWGRbWAcqGHQms8cqGGDbqxcqGHRaWkcqWIMaYsaeaajuaGcaaMc8UaaGPaVpaalaaabaGaeGymaedabaGaeGOmaidaaOWaaiWaaeaacqGGBbWwjuaGdaWcaaqaaiabdchaWjabgMYiHlabdUgaRnaaCaaabeqaaiabikdaYaaacqGHQms8aeaacqGHPms4cqWGRbWAcqGHQms8aaGccqGHsisljuaGdaWcaaqaaiabek7aIbqaaiabeU7aSbaakiabcIcaOiabigdaXiabgkHiTiabdchaWjabcMcaPiabgMYiHlabdUgaRjabgQYiXlabc2faDnaaCaaaleqabaGaeGOmaidaaOGaey4kaSscfa4aaSaaaeaacqaHYoGyaeaacqaH7oaBaaGccqaI0aancqGGOaakcqaIXaqmcqGHsislcqWGWbaCcqGGPaqkcqWGWbaCcqGHPms4cqWGRbWAcqGHQms8daahaaWcbeqaaiabikdaYaaaaOGaay5Eaiaaw2haamaaCaaaleqabaGaeGymaeJaei4la8IaeGOmaidaaaaaaaa@8DF5@
(23)

Remark: From the above figures, we can find that the thresholds are very small, this is because that we are simulating the ratios between two parameters λ, β. In fact, if we fix one of them, and simulate the threshold for another parameter by changing p, here p can go from 0 to 1. For instance, in order to demonstrate this case, we give simulations to verify inequality (19), by letting β = 0.1 unchanged, and consider the threshold for λ, with p changing (see Figure 3).

Figure 3
figure 3

Verification of inequality (19). Simulations show that λ is in accordance with inequality (19) by changing p. And we can see that the smaller p is, the larger the threshold for λ is.

Conclusion

In order to explain the mechanism of the spreading of epidemic, many results were proposed by researchers. These results were often based on that the spreading of epidemics on complex networks via a single mechanism, hence, these results were often not so realistic. Moreover, many results just gave the threshold for the outbreak of epidemics, which can not explicitly answer whether the epidemic will prevail locally or globally.

Because many diseases can spread in different ways, e.g., HIV virus, the avian flu, and so on, in this paper we consider that the spreading of epidemics on complex networks via multiple routes. Under such assumption, we get positive thresholds, neither a constant value nor a zero value, which are more realistic. Furthermore, we answer that the epidemic will prevail globally once some conditions are satisfied.

In Section 3, the results can be obtained by solving a self-consistency equation, and the threshold is the linear combination of λ and β. However, when we consider the case as in Section 4, the results can not be obtained by solving a self-consistency equation, moreover, the threshold is not the linear combination of λ and β, which suggest that the threshold for the multiple routes of spreading epidemic is not just the addition of the respective way of spreading epidemic, sometimes, some complex behavior can happen. Finally, in Section 4, we present numerical simulations to support our analytical results.

References

  1. Kermack WO, McKendrick AG: Contributions to the mathematical theory of epidemics. Proc Roy Soc. 1927, A115: 700-10.1098/rspa.1927.0118.

    Article  ADS  Google Scholar 

  2. Barabási AL, Albert R: Emergence of scaling in random networks. Science. 1999, 286: 509-10.1126/science.286.5439.509.

    Article  ADS  MathSciNet  Google Scholar 

  3. Pastor-Satorras R, Vespignani A: Epidemic dynamics and endemic states in complex networks. Phys Rev E Stat Nonlin Soft Matter Phys. 2001, E63 (Pt 2): 066117-

    Article  Google Scholar 

  4. Pastor-Satorras R, Vespignani A: Epidemic spreading in scale-free networks. Phys Rev Lett. 2001, 86: 3200-10.1103/PhysRevLett.86.3200.

    Article  ADS  Google Scholar 

  5. Pastor-Satorras R, Vespignani A: Infection dynamics in finite size scale-free networks. Phys Rev. 2002, E65: 035108-

    ADS  Google Scholar 

  6. Zhou LJ, Gea T: Behaviors of susceptible-infected epidemics on scale-free networks with identical infectivity. Phys Rev E Stat Nonlin Soft Matter Phys. 2006, E74 (5 Pt 2): 056109-

    Article  ADS  Google Scholar 

  7. Small M, Tse CK: Small world and scale free model of transmission of SARS. Int J Bifurcat and Chaos. 2005, 15: 1745-10.1142/S0218127405012776.

    Article  ADS  Google Scholar 

  8. Small M, Tse CK: Clustering model for tranmsmission of the SARS virus: application to epidemic control and risk assesment. Physica. 2005, A351: 499-

    Article  ADS  Google Scholar 

  9. Masuda N, Konno N: Multi-state epidemic processes on complex networks. J of Theoretical Biology. 2006, 243 (1): 64-74. 10.1016/j.jtbi.2006.06.010.

    Article  MathSciNet  Google Scholar 

  10. Fu XC, Small M, Walker DM, Zhang HF: Epidemic dynamics on scale-free networks with piecewise linear infectivity and immunization. Phys Rev. 2008, E77: 036113-

    ADS  MathSciNet  Google Scholar 

  11. Lajmainovitch A, Yorke JA: A deterministic model for gonorrhea in a nonhomogenous population. Math Biosci. 1976, 28: 221-10.1016/0025-5564(76)90125-5.

    Article  MathSciNet  Google Scholar 

Download references

Acknowledgements

The research was supported jointly by the Foundation of Anhui Education Bureau (Grant KJ2007A003), Natural Science Foundation of Anhui, China (Grant 070416225), a grant from the Health, Welfare and Food Bureau of the Hong Kong SAR Government, and by NSFC Grant 10672146.

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Haifeng Zhang.

Additional information

Authors' contributions

All authors contributed equally to this work.

Authors’ original submitted files for images

Rights and permissions

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and permissions

About this article

Cite this article

Zhang, H., Small, M. & Fu, X. Global behavior of epidemic transmission on heterogeneous networks via two distinct routes. Nonlinear Biomed Phys 2, 2 (2008). https://doi.org/10.1186/1753-4631-2-2

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/1753-4631-2-2

Keywords