# An Epidemic Model of Computer Worms with Time Delay and Variable Infection Rate

Authors: Yu Yao; Qiang Fu; Wei Yang; Ying Wang; Chuan Sheng
Journal: Security and Communication Networks
 (2018)

---

## Abstract

With rapid development of Internet, network security issues become increasingly serious. Temporary patches have been put on the infectious hosts, which may lose efficacy on occasions. This leads to a time delay when vaccinated hosts change to susceptible hosts. On the other hand, the worm infection is usually a nonlinear process. Considering the actual situation, a variable infection rate is introduced to describe the spread process of worms. According to above aspects, we propose a time-delayed worm propagation model with variable infection rate. Then the existence condition and the stability of the positive equilibrium are derived. Due to the existence of time delay, the worm propagation system may be unstable and out of control. Moreover, the thresholdτ0 of Hopf bifurcation is obtained. The worm propagation system is stable if time delay is less than τ0. When time delay is over τ0, the system will be unstable. In addition, numerical experiments have been performed, which can match the conclusions we deduce. The numerical experiments also show that there exists a threshold in the parameter a, which implies that we should choose appropriate infection rate β(t) to constrain worm prevalence. Finally, simulation experiments are carried out to prove the validity of our conclusions.


---

## Body

## 1. Introduction

With the deep application of the Internet, network security plays a more and more important role in recent years. Among the security events, the consequences of large-scale network attacks (such as worm attacks and DOS attacks) are especially serious. Meanwhile, the characteristics of worm attacks are wide infection scale, fast spread speed, and serious harm. Consequently, many experts focus on the spread of Internet worms. Some traditional epidemic models of infectious diseases were used to describe the propagation of Internet worms [1] when the Red Code worms broke out. In order to study the spread of malware among mobile phones, the SIS model [2] is proposed by some researchers. Qing and Wen introduced the Kermack-McKendrick model, which is also called SIR model [3]. Then many mathematical models [4–9] inspired by the SIR model have been employed to constrain the propagation of Internet worms. Some research achievements [10–12] showed that the spread dynamic system of malware would be unstable and bifurcation and chaos would appear. Considering the fact that the intrusion detection system (IDS) may lead to time delay, Yao et al. [13–15] obtained the threshold of time delay when Hopf bifurcation occurred. Pulse quarantine strategy [16] also has been taken to constrain the propagation of worms in network. Due to the effect of different topologies, some experts presented different models [9, 17] to analyse the results.Although most of previous works can offer useful insight into the Internet worm propagations, some of them fail to grasp the detail that has important impact on the worm propagation. Namely, some of the previous models ignore the variation of the infection rate. They usually regard it as a constant that cannot describe the characteristics and dynamics of worm propagation accurately, such as SIS model [18], SIR models [12, 19], SIRS model [20], and SIES model [21]. Moreover, some unconventional models such as delayed models [6, 22] and impulsive models [23, 24] have been proposed. Analogously, these models regard infection rate as a constant as well. In the early stages of worm invasion, the number of infected nodes is small and the linear assumption is still more reasonable. However, as the number of infected nodes increases, the true infection rate will tend to be saturated, and it can be significantly nonlinear. In this case, the linear assumption will overestimate the harmfulness of the worms and lead to great waste of resources.In this paper, a variable infection rate is introduced into the worm propagation. Some experts have suggested that worm infection is a nonlinear process. The majority of previous models mentioned above are based on the bilinear incidence rate assumption, which is a good approximation of the general incidence rate in the case where the proportion of infected computers is small. However, in reality, the density of infected computers may be large [25]. To understand the spreading behaviour of worm propagation better, it is necessary to study epidemic models with general incidence rate. The nonlinear infection rate is used to capture the dynamics of overcrowded infectious networks and high viral loads [26]. Gan et al. [25] show that some nonlinear incidence rates may be conducive to the containment of computer viruses. Feng et al. [27] have proposed the SIRS model with a variable infection rate which plays an important role in the spread of the Internet worm. We consider that the vaccinated hosts (the immunizing hosts) may turn to susceptible hosts (the hosts liable to infection by worms) if the worm variants appear or the patches lose efficacy, and this process may take a period of time. Due to the existence of time delay, vaccinated hosts go through a temporary state (delayed state) after the failure of vaccination before becoming susceptible. In this paper, we try to establish a realistic worm propagation model, motivated by the works [7, 27]. This model can give deep insight into predicting worm spread in networks.The subsequent materials of this work are organized as follows. In Section2, we present the SIQVD model. Section 3 analyses the stability of equilibrium and the threshold of Hopf bifurcation. In Section 4, we carry out the numerical analysis and simulation of our model. Section 5 gives the conclusion and proposes useful strategies.


## 2. Model Formulation

We propose a model of worm propagation to describe the spreading behaviour of Internet worms more realistically in this paper. Susceptible hosts can turn to the infectious state by many factors. Many classical models employ bilinear infection rate described byβSI, where β is determined by the probability of transmission contact betweenS (susceptible hosts) andI (infectious hosts). Previous models usually regard β as a constant. In fact, the worm infection is a nonlinear process so that β should be adjusted to β(t). Infectious hosts can change to vaccinated hosts if there are countermeasures applying to them. The countermeasures include antivirus software, firewall, and patching. Meanwhile, we consider zero-day attacks in this paper. Zero-day attacks spread Internet worms through vulnerabilities of the system or software. Usually, the time of the whole process is not over 24 hours. There are no effective and safe patches when the zero-day attacks appear. So quarantine strategy is proposed to control the worm propagation for the hosts without useful patches. The application of the quarantine strategy relies on the hybrid intrusion detection system (IDS). The hybrid IDS not only can detect unknown worms making up for the lack of misuse detection system but also can avoid the high rate of false positives generated by the anomaly detection system. Thus, some hosts are in quarantined state. And the quarantined hosts can turn to vaccinated hosts by installing patches. Usually some patches are temporary and the temporary patches may lose efficacy if we install operating systems. When the worm variants and unknown worms appear, vaccinated hosts may change to susceptible states. This process can generate a time delay which we called delay states.We assume that all the hosts change over time among five states: susceptible (S), infectious (I), quarantined (Q), vaccinated (V), and delay (D). Let S(t), I(t), Q(t), V(t), and D(t) denote the number of susceptible, infectious, quarantined, vaccinated, and delay hosts, respectively, at time t. We assume that the total number of all the hosts throughout Internet isN. The transition diagram is given in Figure 1.Figure 1
State transition diagram of the SIQVD model.In order to show the parameters clearly, we list some frequent notations of the model in Notations.After above description, we can express the model with the following equations:(1)dSdt=-βtItSt-φSt+δVt-τ,dIdt=βtItSt-γ+αIt,dQdt=αIt-θQt,dVdt=φSt+γIt+θQt-δVt,dDdt=δVt-δVt-τ, where β(t) changes with time t. We regard the infection rate as β(t)=β0f1(I(t)), where f1 is a nonlinear function of I [27]. For being nonlinear, the function f1 is assumed to satisfy the following assumptions [28]:(1)
f1(0)=0.(2)
f1′(I)>0.(3)
f1′′(I)<0.(4)
limt→∞⁡f1I⁡f1I=C<+∞.In other words,f1 is an increasing function that is bounded (by the constant C).From the above discussion, we can express the model by the following differential equations:(2)dSdt=-β0fItSt-φSt+δVt-τ,dIdt=β0fItSt-γ+αIt,dQdt=αIt-θQt,dVdt=φSt+γIt+θQt-δVt,dDdt=δVt-δVt-τ, where f(I(t))=f1(I(t))I(t).


## 3. Stability of Equilibrium and Bifurcation Analysis

Let(3)R0=Nδβ0f′0φ+δα+γ.Theorem 1.
System (2) has a unique positive equilibrium point E∗=(S∗,I∗,Q∗,V∗,D∗) when R0≥1, where(4)S∗=α+γI∗β0fI∗,I∗=I∗,Q∗=αI∗θ,V∗=φ+β0fI∗α+γδβ0fI∗I∗,D∗=δτV∗.Proof.
When system (2) is stable, it satisfies the following equations: (5)-β0fItSt-φSt+δVt-τ=0,β0fItSt-γ+αIt=0,αIt-θQt=0,φSt+γIt+θQt-δVt=0,δVt-δVt-τ=0.
We makeI(t)>0; then we have(6)S∗=α+γI∗β0fI∗,Q∗=αI∗θ,V∗=φ+β0fI∗α+γδβ0fI∗I∗,D∗=δτV∗.
Since the total number of hosts in system (5) is N, we can get the following equation of I:(7)HI=α+γIβ0fI+I+αIθ+φ+β0fIα+γδβ0fII-N=0.
Then we calculate the sign of its derivative as follows:(8)H′I=δ2β0γ+αfI-If′I+δ2β02f2I1+α/θδβ0fI2+δβ02γ+αf2I+δβ0φγ+αfI-If′Iδβ0fI2.
Sincef1′(I)<0, we can get f(I)-If′(I)>0. As a result, H′(I)>0. If there exists a positive root of H(I)=0, H(I) must satisfy H(0)<0. So we can get(9)H0=limI→0⁡HI=δα+γI+δβ0fII+α/θδβ0fII+φ+β0fIα+γI-Nδβ0fIδβ0fI=δ+φα+γδβ0f′01-Nδβ0f′0δ+φα+γ.
From (9), we can conclude that H(0)<0 when R0>1. Hence, there exists a positive equilibrium point if R0>1. The proof is completed.SinceQ(t)=N-S(t)-I(t)-V(t)-D(t), we can simplify system (2) as follows:(10)dStdt=-β0fItSt-φSt+δVt-τ,dItdt=β0fItSt-γ+αIt,dVtdt=φSt+γIt+θN-St-It-Vt-Dt-δVt,dDtdt=δVt-δVt-τ.The Jacobi matrix of system (10) about E∗=(S∗,I∗,Q∗,V∗,D∗) is given by(11)JE∗=β0fI+φβ0f′IS0-δe-λτ-β0fIβ0f′IS-α+γ00000δe-λτ-δθ-φθ-γθ-θ-δ.The characteristic equation of the matrix (11) can be obtained by(12)Pλ+Qλe-λτ=0.The expressions ofP(λ)andQ(λ) are(13)Pλ=λ4+p3λ3+p2λ2+p1λ+p0,Qλ=q2λ2+q1λ+q0,where(14)p3=-β0f′IS+α+γ+β0fI+φ+θ+δ,p2=-β0f′ISθ+δ+φ+β0fIθ+δ+α+γ+α+γθ+δ+φ+θ+δφ+θδ,p1=-β0f′ISθδ+δφ+φθ+β0fIθγ+δγ+αθ+δα+θδ+α+γθδ+δφ+φθ+δθφ,p0=-β0f′ISφ+β0fIα+γ+α+γφδθ,q2=-φδ,q1=β0f′ISφδ-β0fIδγ-α+γδφ-φθδ,q0=--β0f′ISφ+β0fIα+γ+α+γφδθ.Theorem 2.
If the condition is satisfied as(H1), the positive equilibrium E∗=(S∗,I∗,Q∗,V∗,D∗) is locally asymptotically stable without time delay.(15)H1:p3>0,p3p2+q2-p1+q1>0,p1+q1>0.Proof.
Whenτ=0, (12) simplifies to(16)λ4+p3λ3+p2+q2λ2+p1+q1λ+p0+q0=0.
According to Routh-Hurwitz criterion, all the roots of (16) have negative real parts. Hence, we can deduce that the positive equilibrium E∗=(S∗,I∗,Q∗,V∗,D∗) is locally asymptotically stable without time delay. The proof is completed.Assuming thatλ=iw is the root of (12), we substitute it into (12). After separating the real and imaginary parts, we can get the following two equations:(17)w4-p2w2+p0+q0-q2w2cos⁡wτ+q1wsin⁡wτ=0.(18)-p3w3+p1w-q0-q2w2sin⁡wτ+q1wcos⁡wτ=0,Uniting (17) and (18), we can obtain(19)q0-q2w22+q1w2=w4-p2w2+p02+-p3w3+p1w2.It can be written as(20)w8+D3w6+D2w4+D1w2+D0=0, where(21)D3=p32-2p2,D2=2p0+p22-2p3p1-q22,D1=p12-2p2p0+2q0q2-q12,D0=p02-q22.Letz=w2; then (20) can be turned into(22)hz=z4+D3z3+D2z2+D1z+D0.Theorem 3.
Assume that(H1) is satisfied; (a)  D0≥0, Δ≥0, z1>0, or h(z1)>0; (b)  D0≥0, Δ<0, and there is no z∗∈(z1,z2,z3) such that z∗>0 and h(z∗)≤0. Then the positive equilibrium E∗=(S∗,I∗,Q∗,V∗,D∗) of system (1) is absolutely stable. Namely, E∗=(S∗,I∗,Q∗,V∗,D∗) is asymptotically stable for any time delay τ≥0.Assume that the coefficients inh(z) satisfy the following condition:(H3):
(a)D0≥0,Δ≥0,z1>0, or hz1>0; (b) D0≥0,Δ<0, and there is no z∗∈(z1,z2,z3) such that z∗>0 and hz∗≤0.According to previous lemmas, it can be known that (22) has at least a positive root w0, which also means that the characteristic equation (12) has a pair of purely imaginary roots ±iw0.Since the pair of purely imaginary roots±w0 is the roots of (12), we can get the corresponding τk>0 by uniting (17) and (18).(23)τk=1w0arccosq2w2-q0w4-p2w2+p0+p3w3-p1wq1wq2w2-q02+q1w2+2kπw0.Letλ(τ)=v(τ)+iw(τ) be the root of (12). It is satisfied that v(τk)=0,w(τk)=w0.Lemma 4.
Suppose thath′(z0)≠0. If τ=τ0, then ±iw0 is a pair of purely imaginary roots of (12). Moreover, if the conditions in Lemma 3.4(1) in [14] are satisfied, then(24)dRe⁡λτ0dτ>0.This means that there exists at least one eigenvalue with positive real part whenτ>τk. Differentiating on both sides of (12) with respect to τ, we can obtain(25)dλdτ-1=4λ3+3p3λ2+2p2λ+p1+2q2λ+q1e-λτ-q2λ2+q1λ+q0τe-λτq2λ2+q1λ+q0λe-λτ=4λ3+3p3λ2+2p2λ+p1eλτq2λ2+q1λ+q0λ+2q2λ+q1q2λ2+q1λ+q0λ-τλ.According to (17) and (18), we obtain the following:(26)sgndRe⁡λdττ=τk=sgnRedλdτ-1λ=iw0=sgnRe4λ3+3p3λ2+2p2λ+p1eλτq2λ2+q1λ+q0λ+2q2λ+q1q2λ2+q1λ+q0λ-τλλ=iw0=sgnRe-4w03i-3p3w02+2p2iw0+p1cos⁡w0τk+isin⁡w0τk+2q2iw0+q1-q2w02+q1w0i+q0iw0=sgn1Γ4w08+3p32-6p2w06+4p0+2p22-4p3p1-2q22w04+p12-2p2p0+2q0q2-q12w02=sgnw02Γ4w06+3D3w04+2D2w02+D1=sgnw02Γh′w02=sgnh′w02,where Γ=q0w0-q2w032+q1w022.Then it follows hypothesis(H3) and h′(w02)≠0. Therefore(27)dReλdττ=τk>0.According to Routh’s theorem, the root of characteristic equation (12) crosses from left to right on the imaginary axis as τ continuously varies from a value less than τk to one greater than τk. Hence, according to Hopf bifurcation theorem for functional differential equations, the transverse condition holds and the conditions for Hopf bifurcation are satisfied at τ=τk.Theorem 5.
Supposing that the conditionsH1 and H2 are satisfied,(1)
whenτ<τ0, the positive equilibrium  E∗=S∗,I∗,Q∗,V∗,D∗ of system (2) is locally asymptotically stable and it is unstable when τ≥τ0,(2)
when system (2) satisfies H3, the system undergoes a Hopf bifurcation at the positive equilibrium E∗=(S∗,I∗,Q∗,V∗,D∗) when τ=τ0.
This implies that when the time delayτ<τ0, the system will stabilize at its equilibrium point, which is beneficial for us to implement a containment strategy; when the delay τ≥τ0, the system will be unstable and worms cannot be effectively controlled.


## 4. Numerical Simulations and Simulations Experiments

In order to verify the theorems proposed in this paper, we have made the numerical experiments in this section. We select the Slammer worm for experiments. The total number of hostsN is assumed as 400000. Based on the actual situation, the worm’s average scan rate is η=4000 per second. We can calculate the infection rate β=η/232=0.00000093. The susceptible hosts change to vaccinated hosts at rate φ=0.001. The recovered rate of infectious hosts is set as γ=0.002. The quarantine rate α of infectious hosts is 0.2 and the immunity rate θ of quarantined hosts is 0.05. The vaccinated hosts lose immunity at rate δ=0.08. We choose the nonlinear function f1(I(t))=1/(1+aI(t)); then we have f(I(t))=I(t)/(1+aI(t)), where a is the parameter that represents the infection rate sensitivity to the number of infected hosts I(t) [27]. When a is zero, it means that the infection rate is a constant. Then we can get R0=1.819>1. At first, the number of infectious hosts is five and the others’ states are susceptible.Whenτ=10<τ0, we can see the changes of the numbers of four kinds of hosts in Figure 2. From Figure 2, we can find that every kind of hosts will be stable when t=400, which implies that E∗ is locally asymptotically stable. Figure 3 shows the numbers of susceptible, infectious, quarantined, vaccinated, and delayed hosts when τ=60>τ0. In this figure, it can be clearly found that the curves of hosts are fluctuant and it is hard for us to predict the propagation of worms.Figure 2
The worm propagation of the four kinds of hosts’ results withτ<τ0.Figure 3
The worm propagation of the four kinds of hosts’ results withτ>τ0.In order to see the influence of time delay, Figure4 shows the number of infectious hosts in the same coordinate with different time delays τ=5, τ=15, τ=45, and τ=65. Initially, time delay has little effect in the initial stage of worm propagation, which can be obtained by the overlap of the four curves. With the increase of time delay, the curve begins to oscillate. The infecting process gets unstable with time delay passing through the threshold τ0, which meets our conclusions.Figure 4
The number of infected hostsI(t) with the change of τ.Figures5 and 6 show the number of infected hosts in the same condition with a=0.0000001, a=0.000001, and a=0.00001 when τ<τ0 and τ>τ0. From these two figures, we can get the conclusion that the larger a is, the lower peak of the number of infectious hosts is. Therefore, we can choose appropriate a to get proper β(t) to constrain the spread of Internet worms.Figure 5
Relationship betweena and the number of I(t) when τ=5.Figure 6
Relationship betweena and the number of I(t) when τ=85.Figure7 shows the phase portrait of susceptible hosts S(t) and infectious hosts I(t) of system (2) when τ=30<τ0. Moreover, Figure 8 shows the condition when τ=60>τ0. From the figures, we can find that the curve converges to a fixed point, which implies that the system is stable when τ=30<τ0 and the curve radiates to a limit cycle, which implies that the system is unstable when τ=60>τ0. Figures 9 and 10 are the projection of the phase portrait of system (2) in (S,I,V)-space at τ=30 and τ=60. The same conclusion can be obtained by the figures.Figure 7
The phase portrait ofS(t) and I(t) when τ=30.Figure 8
The phase portrait ofS(t) and I(t) when τ=50.Figure 9
The projection of the phase portrait of system (2) in (S,I,V)-space when τ=30.Figure 10
The projection of the phase portrait of system (2) in (S,I,V)-space when τ=50.Figure11 gives the bifurcation diagram of system (2) with the parameter a=0.0000001. It can be easily obtained that the Hopf bifurcation occurs at τ=τ0=37, which is similar to results of theoretical derivation. Figure 12 gives the bifurcation diagram of system (2) with the parameter a=0.00001. The Hopf bifurcation occurs at τ=τ0=76. Comparing the two figures, it is shown that the parameter a has effect on the time of Hopf bifurcation occurrence. As the parameter a increases, the Hopf bifurcation occurs at a later time.Figure 11
Bifurcation diagram of system (2) with a=0.0000001.Figure 12
Bifurcation diagram of system (2) with a=0.00001.In order to simulate the actual behaviour of worm propagation and verify the correctness of the theoretical analysis and numerical simulation, we carry out the discrete-time simulation, which is an expanded version of Zou et al.’s [7] program. The simulation experiment is used to simulate the worm propagation in the real network. There are 400000 hosts in our simulation experiments. At first, we randomly choose five hosts in the network to be infectious hosts and the others’ states are set to be susceptible. In the simulation experiments, the implementation of transition rates of the worm propagation model depends on probability.Figure13 shows the comparisons between numerical and simulation curves of susceptible, infectious, quarantined, and vaccinated hosts when τ=10<τ0, which implies that the simulation curves match the numerical curves very well.Figure 13
Comparisons between numerical and simulation curves of system (2) when τ=10<τ0.
(a) (b) (c) (d)When the value ofτ increases and passes over the threshold value of τ0, namely, τ=60>τ0, numerical and simulation curves of susceptible, infectious, quarantined, and vaccinated hosts can also match very well as Figure 14 shows. We can find that there exists a difference between numerical and simulation curves because of the high precision of numerical and simulation curves because of the high precision of numerical experiment. However, the small difference does not affect the validity of our conclusions.Figure 14
Comparisons between numerical and simulation curves of system (2) when τ=60>τ0.
(a) (b) (c) (d)


## 5. Conclusions

In this paper, we propose a SIQVD model with the variable infection rate based on the consideration of a quarantine strategy. Then we analyse the stability of the positive equilibrium and Hopf bifurcation. The critical time delayτ0 in which Hopf bifurcation appears is obtained. Through the theoretical analysis, the useful conclusions are obtained, which can be verified by the numerical experiments and simulations. The following conclusions can be derived by the current research:(1)
The worm propagation system is stable when time delayτ<τ0. In this condition, we can predict the spread of worms correctly, and the worms can be reduced to a low extent at last.(2)
The worm propagation system is unstable when time delayτ≥τ0, and the system is out of control. Therefore, time delay should be controlled in a proper range: τ<τ0.(3)
When parametera increases an order of magnitude, the infection rate reduces and the peak of the number of infected hosts will decrease very obviously. Meanwhile, R0 will be reduced by the decrease of f′(0). Hence, there exists a threshold for a, and we can choose proper infection rate by adjusting the value of a to control the prevalence of worms.The worm propagation model can be used for Internet worms, such as Code Red worms, Slammer worms, and Witty worms. It can predict the spreading behaviour of Internet worms more realistically. In our future work, we will focus more on the network structure and study it further.



---

*Source: 9756982-2018-03-06.xml*