# Dynamical Analysis of a Computer Virus Propagation Model with Delay and Infectivity in Latent Period

Authors: Zizhen Zhang; Dianjie Bi
Journal: Discrete Dynamics in Nature and Society
 (2016)

---

## Abstract

A delayed SLB computer virus propagation model with infectivity in latent period is proposed in this paper. We establish sufficient conditions for local stability of the positive equilibrium and existence of Hopf bifurcation by analyzing distribution of the roots of the associated characteristic equation and applying the Hopf bifurcation theorem. Furthermore, properties of the Hopf bifurcation are determined by using the normal form theory and the center manifold theorem. Finally, numerical simulations supporting the theoretical analysis are carried out.


---

## Body

## 1. Introduction

Due to their striking features such as destruction, polymorphism, and unpredictability [1–3], computer viruses, ranging from conventional viruses to network worms and other malicious codes such as Trojans and spyware, have addressed a serious threat to confidentiality, integrity, and availability of computer resources on the Internet. They have also come to be one major threat to our work and daily life.In recent years, people have paid attention to the necessity of monitoring computer viruses in the Internet and proposed many mathematical models to depict the spread of computer viruses in the Internet based on the classic epidemic models [4–9]. In [4], Li et al. proposed a worm propagation model with removable storage device and studied stability of the model. Mishra and Keshri proposed a SEIRS-V model to describe the dynamics of worm propagation in wireless sensor networks [6]. However, most of the above models assume that only breaking-out computers have infectivity. This is not consistent with reality. Because computer viruses have a feature of latency, an infected computer which is in latency can also infect other computers through file copying or file downloading. Based on this consideration, in [8], Yang et al. proposed the following computer virus model with graded cure rate:(1)dStdt=δ-βStLt+Bt+γ1Lt+γ2Bt-δSt,dLtdt=βStLt+Bt-γ1Lt-αLt-δLt,dBtdt=αLt-γ2Bt-δBt,where S(t) is the percentage of uninfected computers in the Internet, at time t. L(t) and B(t) are the percentages of latent and breaking-out computers in the Internet, at time t, respectively. δ is the connected rate of external computers to the Internet and it is also the disconnected rate of internal computers from the Internet. β is the infected rate of the uninfected computers. α is the breaking-out rate of the latent computers. γ1 and γ2 are the cured rates of the latent and breaking-out computers, respectively. Obviously, it is assumed that both seizing and latent computers have infectivity.However, system (1) ignores time delay due to the period that a computer user uses to cure the latent and the breaking-out computers by using antivirus software or reinstalling system. As is known, time delay can make a dynamical system lose its stability and then lead to the occurrence of a Hopf bifurcation. Hopf bifurcation of dynamical systems with delay has been investigated by many authors [10–16]. In [10], Zhuang and Zhu analyzed the existence of Hopf bifurcation for an improved HIV model with time delay and cure rate by regarding the time lag from infection of cells to the cells becoming actively infected as a bifurcation parameter. In [11], Liu investigated the existence and properties of the Hopf bifurcation in a delayed SEIQRS model for the transmission of malicious objects in computer network. In [14], Bianca et al. studied the bifurcation of a delayed-energy-based model of capital accumulation. Motivated by the work above, we introduce the time delay due to the period that a computer user uses to cure the latent and breaking-out computers into system (1) and get the following delayed computer virus model:(2)dStdt=δ-βStLt+Bt+γ1Lt-τ+γ2Bt-τ-δSt,dLtdt=βStLt+Bt-γ1Lt-αLt-δLt,dBtdt=αLt-γ2Bt-δBt,where τ is the delay due to the period that a computer user uses to cure the latent and breaking-out computers by using antivirus software or reinstalling system.The remainder of this paper is organized as follows. Local stability and existence of Hopf bifurcation for system (2) are analyzed in Section 2. Properties of the Hopf bifurcation are investigated in Section 3 by means of the center manifold theorem and the normal form theory. A numerical example is presented in Section 4 to illustrate the theoretical results. This work is summarized in Section 5.


## 2. Existence of Hopf Bifurcation

According to the analysis in [8], we can conclude that if the basic reproduction number R0=β(α+γ2+δ)/(γ2+δ)(α+γ1+δ)>1, then system (2) has a unique positive equilibrium E+(S∗,L∗,B∗), where(3)S∗=α+γ1+δL∗βL∗+B∗,L∗=γ2+δ1-1/R0α+γ2+δ,B∗=α1-1/R0α+γ2+δ.For the linearized system of system (2) at the positive equilibrium E+(S∗,L∗,B∗), the corresponding characteristic equation is(4)λ3+A2λ2+A1λ+A0+B1λ+B0e-λτ=0, where(5)A0=a11a23a32-a22a33+a21a12a33-a13a32,A1=a11a22+a33+a22a33-a12a21-a23a32,A2=-a11+a22+a33,B0=a21a33b12-a32b13,B1=-a21b12,a11=-βL∗+B∗-δ,a12=-βS∗,a13=-βS∗,a21=βL∗+B∗,a22=βS∗-γ1+α+δL∗,a23=βS∗,a32=α,a33=-γ2+δ,b12=γ1,b13=γ2.When τ=0, (4) becomes(6)λ3+A2∗λ2+A1∗λ+A0∗=0,where(7)A0∗=B0+C0,A1∗=A1+B1,A2∗=A2.Obviously, if condition(H1), A2∗>0, A2∗A1∗>A0∗>0, holds, then system (2) is locally asymptotically stable when τ=0. Next, we choose the time delay τ as the parameter to consider the local stability of the positive equilibrium E+(S∗,L∗,B∗) and the Hopf bifurcation of system (2).Forτ>0, let λ=iω(ω>0) be the root of (4); then(8)B1ωcosτω-B0sinτω=ω3-A1ω,B1ωsinτω+B0cosτω=A2ω2-A0,from which we can get the following sixth-degree equation for ω:(9)ω6+c2ω4+c1ω2+c0=0,where(10)c2=A22-2A1,c1=A12-2A0A2-B12,c0=A02-B02.Let ω2=v; then (9) becomes(11)v3+c2v2+c1v+c0=0.From the analysis above, we know that if the coefficients of system (2) are given, one can easily get the coefficients of (11) by MATLAB software package. Correspondingly, the roots of (11) can be obtained. Thus, we make the following assumption in order to give the main results in this paper.( H 2 ) Equation (11) has at least one positive root.If condition(H2) holds, then there exists one positive root v0 of (11) such that (4) has a pair of purely imaginary roots ±iω0=±iv0. For ω0, the corresponding critical value of delay is(12)τ0=1ω0arccosB1ω04+A2B0-A1B1ω02-A0B0B02+B12ω02.And then, we will verify the transversality condition. From (4), we can obtain(13)dλdτ-1=-3λ2+2A2λ+A1λλ3+A2λ2+A1λ+A0+B1λB1λ+B0-τλ.This gives(14)Redλdττ=τ0-1=f′v∗B02+B12ω02,where f(v)=v3+c2v2+c1v+c0 and v∗=ω02.It is clear that if condition(H3)f′(v∗)≠0 holds, then Re[dλ/dτ]τ=τ0-1≠0. According to the Hopf bifurcation theorem in [15], we have the following results for system (2).Theorem 1.
If conditions(H1)–(H3) hold, then the positive equilibrium E+(S∗,L∗,B∗) of system (2) is locally asymptotically stable for τ∈[0,τ0) and system (2) undergoes a Hopf bifurcation at the positive equilibrium E+(S∗,L∗,B∗) when τ=τ0 and a family of periodic solutions bifurcate from the positive equilibrium E+(S∗,L∗,B∗) near τ=τ0.


## 3. Properties of the Hopf Bifurcation

In the previous section, we have obtained the conditions under which a Hopf bifurcation occurs and a family of periodic solutions bifurcate from the positive equilibriumE+(S∗,L∗,B∗) of system (2) when the delay passes through the critical value τ0. In this section, we will derive the explicit formulae determining the direction of the Hopf bifurcation and the stability of bifurcating periodic solution of system (2) at τ=τ0. The approach employed here is the normal form method and center manifold theorem introduced by Hassard et al. in [17].Letτ=τ0+μ, u1=S(τt), u2=L(τt), u3=B(τt), μ∈R, Lμ:C→R2, and F:R×C→R3, so that system (2) is transformed into a functional differential equation in C=C([-1,0],R3) as(15)u˙t=Lμut+Fμ,ut,where(16)Lμφ=τ0+μA∗φ0+B∗φ-1,F1φ,μ=-βτ0+μφ10φ20+φ10φ30,F2φ,μ=βτ0+μφ10φ20+φ10φ30,F3φ,μ=0with(17)A∗=a11a12a13a22a22a230a32a33,B∗=0b12b13000000.Then,Lμ is a bounded linear operator in C=C([-1,0],R3). By the Riesz representation theorem, there exists a function η(θ,μ) of bounded variation for θ∈[-1,0] such that(18)Lμϕ=∫-10dηθ,μϕθ,ϕ∈C-1,0,R3.In fact, we choose(19)ηθ,μ=τ0+μA∗δθ+B∗δθ+1,where δ(θ) is the Dirac delta function.Forϕ∈C([-1,0],R3), we define(20)Aμϕ=dϕθdθ,-1≤θ<0,∫-10dηθ,μϕθ,θ=0,Rμϕ=0,-1≤θ<0,Fμ,ϕ,θ=0.Then system (15) is equivalent to the following operator equation:(21)u˙t=Aμut+Rμut.Next, we define the adjoint operator A∗ of A,(22)A∗φ=-dφsds,0<s≤1,∫-10dηTs,0φ-s,s=0,and a bilinear inner product,(23)φs,ϕθ=φ¯0ϕ0-∫θ=-10∫ξ=0θφ¯ξ-θdηθϕξdξ,where η(θ)=η(θ,0).By the discussion in the previous section, we know that±iω0τ0 are common eigenvalues of A(0) and A∗(0). We need to compute the eigenvectors of A(0) and A∗(0) corresponding to +iω0τ0 and -iω0τ0, respectively. Let q(θ)=(1,q2,q3)Teiω0τ0θ be the eigenvector of A(0) corresponding to +iω0τ0 and let q∗(θ)=D(1,q2∗,q3∗)eiω0τ0s be the eigenvector of A∗(0) corresponding to -iω0τ0. Then we have(24)A0qθ=iω0τ0qθ,A∗0q∗Tθ=-iω0τ0q∗Tθ.From the definitions of A(0) and A∗(0), we have(25)A0qθ=dqθdθ,A∗0q∗Ts=-dq∗Tsds.Thus,(26)qθ=q0eiω0τ0θ,q∗s=q0eiω0τ0s.In addition,(27)∫-10dηθqθ=A∗q0+B∗q-1=A0q0=iω0τ0q0.Thus, we can get(28)q2=iω0-a11-a13+b13e-iω0τ0q3a12+b12e-iω0τ0,q3=a32iω0-a11Δ1+Δ2,where(29)Δ1=a32a13+b13e-iω0τ0,Δ2=iω0-a33a12+b12e-iω0τ0.On the other hand, (30)∫-10dηθq∗-s=A1Tq∗T0+B1Tq∗T-1=A∗0q∗T0=-iω0τ0q∗T0.Thus, we get(31)q2∗=-iω0+a11a21,q3∗=iω0+a11iω0+a22a21a32+a12+b12eiω0τ0a32.From (23), we have(32)q∗s,qθ=D¯q¯∗0q0-∫-10∫0θD¯q¯∗0e-iω0τ0ξ-θdηθq0eiω0τ0ξdξ=D¯1+q2q¯2∗+q3q¯3∗-∫-10q¯∗0θeiω0τ0θdηθq0=D¯1+q2q¯2∗+q3q¯3∗+q¯∗0B∗e-iω0τ0q0=D¯1+q2q¯2∗+q3q¯3∗+τ0e-iω0τ0b12q2+b13q3.Thus, we can obtain(33)D¯=1+q2q¯2∗+q3q¯3∗+τ0e-iω0τ0b12q2+b13q3-1,such that 〈q∗,q〉=1.On the other hand, since〈φ,Aϕ〉=〈A∗φ,ϕ〉, we have(34)-iω0τ0q∗,q¯=q∗,Aq¯=A∗q∗,q¯=-iω0τ0q∗,q¯=iω0τ0q∗,q¯.Obviously, 〈q∗,q¯〉=0.In the remainder of this section, we use the same notations as those used by Xu and He [15] and we first compute the coordinates to describe the center manifold C0 at μ=0. Let xt be the solution of (21) when μ=0.Define(35)zt=q∗,xt,Wt,θ=xtθ-2Reztqθ,on the center manifold C0, and we have(36)Wt,θ=Wzt,z¯t,θ,where(37)Wzt,z¯t,θ=Wz,z¯=W02z22+W11zz¯+W02z¯22+⋯,and z and z¯ are local coordinates for center manifold C0 in the direction of q∗ and q¯∗. Note that W is real if xt is real; we only deal with real solutions. It is easy to see that(38)z˙t=q∗,x˙t=q∗,A0xt+q∗,R0xt=A∗0q∗,xt+q∗,R0xt=A∗0q∗,xt+q¯∗0R0xt-∫-10∫0θq¯∗ξ-θdηθA0R0xtξdξ=iω0τ0zt+q¯∗0f0,xtθ≔iω0τ0zt+q¯∗0f0zt,z¯t.We rewrite this equation as(39)z˙t=iω0τ0zt+gz,z¯, where(40)gz,z¯=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋯.So, we can get(41)gz,z¯=q¯∗0f0z,z¯=D¯1,q¯2∗,q¯3∗f10,xt,f20,xt,0T,where(42)f10,xt=-βτ0φ10φ20+φ10φ30,f20,xt=βτ0φ10φ20+φ10φ30.Since(43)xt=xt+θ=Wz,z¯,θ+zqθ+z¯q¯θ,qθ=1,q2,q3Teiω0τ0θ,we have(44)xt=x1t+θx2t+θ0+W1t+θW2t+θ0+z1q2q3eiω0τ0θ+z¯1q¯2q¯3e-iω0τ0θ,φ10=z+z¯+W2010z22+W1110zz¯+W021z¯22+⋯,φ20=zq2+z¯q¯2+W2020z22+W1120zz¯+W022z¯22+⋯,φ30=zq3+z¯q¯3+W2030z22+W1130zz¯+W023z¯22+⋯.From (40) and (41), we have(45)gz,z¯=D¯1,q¯2∗,q¯3∗×K11z2+K12zz¯+K13z¯2+K14z2z¯K21z2+K22zz¯+K23z¯2+K24z2z¯0+⋯,where(46)K11=-βτ0q2+q3,K12=-βτ0q2+q¯2+q3+q¯3,K13=-βτ0q¯2+q¯3,K14=-βτ0W1110q2+12W2010q¯2+W1120+12W2020+W1110q3+12W2010q¯3+W1130+12W2030,K21=βτ0q2+q3,K22=βτ0q2+q¯2+q3+q¯3,K23=βτ0q¯2+q¯3,K24=βτ0W1110q2+12W2010q¯2+W1120+12W2020+W1110q3+12W2010q¯3+W1130+12W2030.Thus, we have(47)gz,z¯=D¯K11+q¯2∗K21z2+K12+q¯2∗K22zz¯+K13+q¯2∗K23z¯2+K14+q¯2∗K24z2z¯+⋯.Comparing the coefficients in (47) with those in (40), we can get(48)g20=2D¯K11+q¯2∗K21,g11=D¯K12+q¯2∗K22,g02=2D¯K13+q¯2∗K23,g21=2D¯K14+q¯2∗K24.In order to get the expression ofg21, we need to compute W11(θ) and W20(θ). From (21) and (39), we have(49)W′=AW-2Req¯∗0fqθ,-1≤θ<0,AW-2Req¯∗0fqθ+f,θ=0.Define (49) as(50)W′=AW+Hz,z¯,θ,where(51)Hz,z¯,θ=H20θz22+H11θzz¯2+H02θz¯22+⋯.From (37), (39), (50), and (51), we obtain(52)2iω0τ0-AW20θ=H20θ,(53)AW11θ=-H11θ.From (40) and (49), we have(54)Hz,z¯,θ=-2Req¯∗0f0qθ=-2Regz,z¯qθ=-gz,z¯qθ-g¯z,z¯q¯θ=-g02z22+g11zz¯+g02z¯22+g21z2z¯2qθ-g¯02z¯22+g¯11zz¯+g¯02z22+g¯21zz¯22q¯θ-⋯.Comparing the coefficients of (51) and (54), we get(55)H20θ=-g20qθ-g¯02q¯θ,(56)H11θ=-g11qθ-g¯11q¯θ,when θ∈[-1,0).From (52) and (55) and the definition of A(0), we have(57)W˙20θ=2iω0τ0W20θ+g20qθ+g¯02q¯θ.According to q(θ)=q(0)eiω0τ0, we get(58)W20θ=ig20ω0τ0q0eiω0τ0+ig¯023ω0τ0q¯0e-iω0τ0+E1e2iω0τ0,and by a similar method, we can get(59)W11θ=-ig11ω0τ0q0eiω0τ0+ig¯11ω0τ0q¯0e-iω0τ0+E2,where E1 and E2 are both constant vectors.In what follows, we will find outE1 and E2. It follows from the definition of A(0) and (55) and (56) that(60)∫-10dηθW20θ=2iω0τ0W200-H200,(61)∫-10dηθW11θ=-H110,where η(θ)=η(0,θ).From (52) and (53), we get(62)H200=-g20qθ-g¯02q¯0+K11,K21,0T,(63)H110=-g11qθ-g¯11q¯0+K12,K22,0T.Notice that(64)iω0τ0I-∫-10dηθeiω0τ0q0=0,-iω0τ0I-∫-10dηθe-iω0τ0q¯0=0,and substituting (58) and (62) into (60), we obtain(65)2iω0τ0I-∫-10dηθe2iω0τ0θE1=K11K210.That is,(66)2iω0τ0I-A∗-B∗e-2iω0τ0E1=K11K210.Thus, we have(67)E1=2τ0a11′a12′a13′-a21a22′-a230-a32a33′-1×K11K210,where(68)a11′=2iω0-a11,a12′=-a12-b12e-2iω0τ0,a13′=-a13-b13e-2iω0τ0,a22′=2iω0-a22,a33′=2iω0-a33.Similarly, we have(69)∫-10dηθE2=K12K220.Then,(70)A∗+B∗E2=K12K220.Therefore, we can get(71)E2=1τ0a11b12′b13′a21a22a230a32a33-1×K12K220,where(72)b12′=a12+b12,b13′=a13+b13.Thus, we can compute the following coefficients:(73)C10=i2ω0τ0g11g20-2g112-g0223+g212,μ2=-ReC10Reλ′τ0,β2=2ReC10,T2=-ImC10+μ2Imλ′τ0ω0τ0,which determine the properties of the bifurcating periodic solutions in the center manifold at the critical value τ0 and we obtain the following results.Theorem 2.
For system (2), (i)
the sign ofμ2 determines the direction of the Hopf bifurcation: if μ2>0(μ2<0), then the Hopf bifurcation is supercritical (subcritical);(ii)
the sign ofβ2 determines the stability of the bifurcating periodic solutions: if β2<0(β2>0), then the bifurcating periodic solutions are stable (unstable);(iii)
the sign ofT2 determines the period of the bifurcating periodic solutions: if T2>0(T2<0), then the bifurcating periodic solutions increase (decrease).


## 4. Numerical Simulation

In this section, we try to present some numerical simulations for system (2) to validate the previous main results. By extracting some values from [8] and considering the conditions for the existence of the Hopf bifurcation, we choose a set of parameters as follows: α=0.3, β=0.35, γ1=0.1, γ2=0.3, and δ=0.05. Then, we can get a specific case of system (2).By a simple computation, we getR0=1.4444>1 and the unique positive equilibrium E+(0.6924, 0.1657, 0.1420) of the system. Further, from (11) and (12), we obtain ω0=0.2766 and τ0=22.9208. According to Theorem 1, we can conclude that if τ∈[0,6.0897), then the positive equilibrium E+(0.6924, 0.1657, 0.1420) is locally asymptotically stable. As can be seen from Figures 1 and 2, when we choose τ=18.3656<τ0, the positive equilibrium E+(0.6924, 0.1657, 0.1420) is locally asymptotically stable. Once the value of the delay τ passes through the critical value τ0, as is shown in Figures 3 and 4, when we choose τ=32.8965>τ0, a Hopf bifurcation occurs and a family of periodic solutions bifurcate from the positive equilibrium E+(0.6924, 0.1657, 0.1420).Figure 1
E + is locally asymptotically stable when τ=18.3656<τ0 with initial value “0.92,0.155,0.112.”Figure 2
The phase plot of the statesS∗, L∗, and B∗ when τ=18.3656<τ0 with initial value “0.92,0.155,0.112.”Figure 3
E + loses stability and a Hopf bifurcation occurs when τ=32.8965>τ0 with initial value “0.92,0.155,0.112.”Figure 4
The phase plot of the statesS∗, L∗, and B∗ when τ=32.8965>τ0 with initial value “0.92,0.155,0.112.”In addition, we also find that onset of the Hopf bifurcation can be delayed if the values of the parametersγ1 and γ2 increase from the numerical simulation. And this can be realized by means of strengthening the immunization of the new computers connected to the network. Thus, we can conclude that the managers of the network should strengthen the immunization of the new computers connected to the network so as to predict and control the propagation of the computer virus in the network easily.


## 5. Conclusion

A new epidemic model of computer virus with time delay has been proposed in this paper. Compared with the model considered in the literature [8], we mainly investigate the effect of the delay due to the period that a computer user uses to cure the latent and breaking-out computers by using antivirus software or reinstalling system on system (2) by analyzing the corresponding characteristic equation and regarding the time delay as a bifurcation parameter.It is found that when the value of the time delay is below the critical valueτ0 the model is asymptotically stable. In this case, the prevalence of the viruses in the Internet can be easily controlled. However, when the value of the time delay is above τ0, the model will lose its stability and a Hopf bifurcation occurs, which indicates that computers of the three classes in system (2) may coexist in an oscillatory mode under some conditions. This phenomenon is not welcome in the Internet and it is difficult to control and prevent the viruses in the Internet in this case. Hence, we should control the occurrence of the Hopf bifurcation by combining some bifurcation control strategies and other relative features of virus prevalence, in order to make the propagation of computer viruses possible to be predicted and controlled easily. Specially, with the help of the normal form theory and the center manifold theorem, explicit formulae determining the properties of the Hopf bifurcation have been derived. Finally, numerical examples are given to illustrate the main results obtained in the paper.



---

*Source: 3067872-2016-12-22.xml*