# Stability and Hopf Bifurcation for a Delayed Computer Virus Model with Antidote in Vulnerable System

Authors: Zizhen Zhang; Yougang Wang; Massimiliano Ferrara
Journal: Journal of Control Science and Engineering
 (2017)

---

## Abstract

A delayed computer virus model with antidote in vulnerable system is investigated. Local stability of the endemic equilibrium and existence of Hopf bifurcation are discussed by analyzing the associated characteristic equation. Further, direction of the Hopf bifurcation and stability of the bifurcating periodic solutions are investigated by using the normal form theory and the center manifold theorem. Finally, numerical simulations are presented to show consistency with the obtained results.


---

## Body

## 1. Introduction

Applications based on computer networks are becoming more and more popular in our daily life. While bringing convenience to us, computer networks are exposed to various threats [1, 2]. Therefore, it is urgent to explore the spreading law of computer viruses in networks. To this end, many dynamical models describing propagation of computer viruses have been established by scholars at home and abroad. Particularly the classic epidemic models, such as SIRS [3–5] model, SEIRS model [6], and SEIQRS model [7, 8], are used to investigate the spreading law of computer viruses due to the common feature between the computer virus and the biological virus. Some computer virus models with infectivity in both seizing and latent computers have been also proposed by Yang et al. [9–13].Recently, Khanh and Huy [14] proposed the following computer virus model with different antidote rates of nodes in vulnerable system considering the immunizations ways and the vulnerabilities of the operating system:(1)dStdt=A-βStLt+ϖIt-εRt-μSt,dLtdt=βStLt-δ+γ+μLt,dItdt=γLt-α+μ+ϖIt,dRtdt=δLt+αIt-μ-εRt,where St, Lt, It, and Rt are the sizes of susceptible, exposed, infectious, and recovered nodes at time t, respectively; A is the constant recruitment of the susceptible nodes; μ is the same rate at which every node in the states St, Lt, It, and Rt disconnects from the network; ε is the constant rate at which every susceptible node acquires temporary immunity due to antidote and Khanh and Huy [14] assume that ε<μ taking system vulnerability into account; α, β, δ, and ϖ are the other state transition rates of system (1). Khanh and Huy [14] studied stability of system (1) and suggested some effective strategies for eliminating viruses.It is well known that time delays of one type or another have been incorporated into computer virus models due to latent period of virus, temporary immunization period of nodes, or other reasons. Computer virus models with time delay have been investigated extensively in recent years [15–18]. Time delays can play a complicated role in the dynamics of dynamical systems, especially that they can cause Hopf bifurcation phenomenon of the models. In reality, the occurrence of Hopf bifurcation means that the state of computer virus spreading changes from an equilibrium point to a limit cycle, which is not welcomed in networks. Motivated by the work above and considering that it needs a period to reinstall system for the infected nodes in the network, we propose the following system with time delay:(2)dStdt=A-βStLt+ϖIt-τ-εRt-μSt,dLtdt=βStLt-δ+γ+μLt,dItdt=γLt-α+μIt-ϖIt-τ,dRtdt=δLt+αIt-μ-εRt,where τ is the time delay due to the period that the infected nodes use to reinstall system.The rest of the paper is organized as follows. Section2 obtains the basic reproduction number of system (2) and discusses local stability of the endemic equilibrium and existence of Hopf bifurcation; Section 3 determines direction of the Hopf bifurcation and stability of the bifurcated periodic solutions; Section 4 covers numerical analysis and simulations. Finally, Section 5 summarizes this work.


## 2. Existence of Hopf Bifurcation

By a direct computation, we get that ifR0=βA/(μδ+γ+μ)>1, then system (2) has a unique endemic equilibrium P∗S∗,L∗,I∗,R∗, where S∗=(δ+γ+μ)/β, L∗=μ-εα+μ+ϖG/(βF), I∗=γμ-εG/(βF), and αδ+αγ+δμ+δϖG/(βF) with(3)F=μαδ+γ+μ-ε+δμ+ϖ+μ-εγ+μ+ϖ,G=βA-μδ+γ+μ,and R0 is the basic reproduction number of system (2).Based on the leading matrix of system (2) at the endemic equilibrium P∗S∗,L∗,I∗,R∗, the characteristic equation can be obtained as the following form:(4)λ-α11-α12-β13e-λτ-α14-α21λ-α22000-α32λ-α33-β33e-λτ00-α42-α43λ-α44=0,which derives that(5)Mλ+Nλe-λτ=0,where(6)Mλ=λ4+α3λ3+α2λ2+α1λ+α0,(7)Nλ=β3λ3+β2λ2+β1λ+β0,with(8)α0=α11α22-α12α21α33α44+α33α42-α32α43α14α21,α1=α12α21-α11α22α33+α44-α11+α22α33α44-α14α21α42,α2=α11α22-α12α21+α33α44+α11+α22α33+α44,α3=-α11+α22+α33+α44,β0=α11α22-α12α21α44β33+α21α14α42β33+α32α44β13,β1=β33α12α21-α11α22-α44β33α11+α22-α21α32β13,β2=β33α11+α22+α44,β3=-β33,α11=-βL∗-μ,α12=-βS∗,α14=ϖ,α21=βL∗,α22=βS∗-δ+γ+μ,α32=γ,α33=-α+μ,α42=δ,α43=α,α44=-μ-ε,β13=ϖ,β33=-ϖ.When τ=0, (5) takes the following form:(9)λ4+α3+β3λ3+α2+β2λ2+α1+β1λ+α0+β0=0.According to the expressions of αi and βi (i=0,1,2,3), we can get(10)α0+β0=μμ-εα+μ+ϖαγ+μ+δ-ε+δμ+ϖ+μ-εγ+μ+ϖQ,α1+β1=μ-εα+μ+ϖμ+2μμ-ε+ϖ+2μ-εα+γ+ϖ+αδ+γ+δϖQ,α2+β2=α2μ-ε+μ-ε2μ+ϖ+μμ+ϖ+μ-εα+μ+ϖα+δ-ε+γ+3μ+ϖQ,α3+β3=α+ϖ+3μ-ε+μ-εα+μ+ϖQ,with Q=G/F and F, G are specified by (6) and (7), respectively. According to the assumption ε<μ considering system vulnerability, we know that α0+β0>0, α1+β1>0, α2+β2>0, and α3+β3>0. In addition, one can be conclude that(11)α1+β1α2+β2α3+β3-α1+β12-α0+β0α3+β32>0based on the analysis by Khanh and Huy in [14] when R0>1. Therefore, we can conclude that the endemic equilibrium P∗S∗,L∗,I∗,R∗ is locally stable for R0>1 according to the Routh-Hurwitz criterion.The time delay is always positive in our physical system. It follows that if instability occurs for a particular value of the delayτ>0, a characteristic root of (5) must intersect the imaginary axis according to Corollary 2.4 in [19]. Assume that (5) has a purely imaginary root iωω>0. Then, the following identity must be true(12)β1ω-β3ω3sin⁡τω+β0-β2ω2cos⁡τω=α2ω2-ω4-α0,β1ω-β3ω3cos⁡τω-β0-β2ω2sin⁡τω=α3ω3-α1ω,which yields the following equation:(13)ω8+p3ω6+p2ω4+p1ω2+p0=0,where(14)p0=α02-β02,p1=α12-2α0α2-β12+2β0β2,p2=α22+2α0-2α1α3+2β1β3-β22,p3=α32-β32-2α2.Denote ω2=v; (13) becomes(15)v4+p3v3+p2v2+p1v+p0=0.Define(16)fv=v4+p3v3+p2v2+p1v+p0.Thus,(17)f′v=4v3+3p3v2+2p2v+p1.Set(18)4v3+3p3v2+2p2v+p1=0.Let y=v+3p3/4. Then, (18) becomes(19)y3+r1y+s1=0,where(20)r1=p22-316p32,s1=p3332-p2p38+p1.Denote(21)D=s122+r133,σ=-1+3i2,y1=-s12+D3+-s12-D3,y2=-s12+D3σ+-s12-D3σ2,y3=-s12+D3σ2+-s12-D3σ,vi=yi-3p34,i=1,2,3.Discussion about the distribution of the roots of (15) is similar to that in [20]. Thus, we have the following lemma.Lemma 1.
For (15), one has the following:(i)
Ifp0<0, (15) has at least one positive root.(ii)
Ifp0≥0 and D≥0, (15) has positive roots if and only if v1>0 and fv1<0.(iii)
Ifp0≥0 and D<0, (15) has positive roots if and only if there exists at least one v∗∈v1,v2,v3 such that v∗>0 and fv∗≤0.In what follows, we assume thatH1and the coefficients in fv satisfy one of the following conditions in (a)–(c).(a)
p0<0.(b)
p0≥0, D≥0, v1>0, and fv1<0.(c)
p0≥0 and D<0, and there exists at least one v∗∈v1,v2,v3 such that v∗>0 and fv∗≤0.Suppose that the conditionH1 holds; then (15) has at least one positive root v0 such that (13) has a positive root ω0=v0. Further, we can get(22)τ0=1ω0arccos⁡β2-α3β3ω06+α1β3+α3β1-α2β2ω04+α0β2+α2β0-α1β1ω02-α0β0β1ω0-β3ω032+β0-β2ω022.Differentiating both sides of (5) with respect to τ yields(23)dλdτ-1=-4λ3+3α3λ2+2α2λ+α1λλ4+α3λ3+α2λ2+α1λ+α0+3β3λ2+2β2λ+β1λβ3λ3+β2λ2+β1λ+β0-τλ.Further, we have(24)Redλdττ=τ0-1=f′v∗∗β1ω0-β3ω032+β0-β2ω022,where fv=v4+p3v3+p2v2+p1v+p0 and v∗∗=ω02.Obviously, if the conditionH2, f′ω02≠0, holds, then Redλ/dττ=τ0-1≠0. Therefore, the transversality condition is satisfied if H2, f′ω02≠0, holds. According to the Hopf bifurcation theorem in [21], we can obtain the following results.Theorem 2.
Suppose that conditionsH1 and H2 hold for system (2). The endemic equilibrium P∗S∗,L∗,I∗,R∗ is locally asymptotically stable when τ∈0,τ0 and a Hopf bifurcation occurs at the endemic equilibrium P∗S∗,L∗,I∗,R∗ when τ=τ0.


## 3. Properties of the Hopf Bifurcation

Letτ=τ0+μ, μ∈R. By the transformation, x1t=St-S∗, x2t=Lt-L∗, x3t=It-I∗, x4t=Rt-R∗, t→t/τ, and system (2) is equivalent to the following in C=C-1,0,R4:(25)x˙t=Lμxt+Fμ,xt,where xt=x1t,x2t,x3t,x4tT=S,L,I,RT∈R4, xtθ=xt+θ∈C, and Lμ, Fμ,xt are defined as follows, respectively:(26)Lμϕ=τ0+μMϕ0+Nϕ-1,Fμ,ϕ=τ0+μ-βϕ10ϕ20βϕ10ϕ2000,with(27)M=α11α120α14α21α22000α32α3300α42α43α44,N=00β130000000β3300000.Thus, there exists a4×4 matrix function whose components are functions ηθ of the bounded variation in θ∈-1,0 such that(28)Lμϕ=∫-10dηθ,μϕθ,ϕ∈C.In fact, we choose(29)ηθ,μ=τ0+μMδθ+Nδθ+1.For ϕ∈C-1,0,R4, define(30)
A
μ
ϕ
=
d
ϕ
θ
d
θ
,
-
1
≤
θ
<
0
,
∫
-
1
0
d
η
θ
,
μ
ϕ
θ
,
θ
=
0
,
R
μ
ϕ
=
0
,
-
1
≤
θ
<
0
,
F
μ
,
ϕ
,
θ
=
0
.Then system (25) is equivalent to(31)x˙t=Aμut+Rμut.For φ∈C10,1,R4∗, we define the adjoint operator A∗ of A0(32)A∗φ=-dφsds,0<s≤1,∫-10dηTs,0φ-s,s=0,and the following bilinear inner product(33)φs,ϕθ=φ¯0ϕ0-∫θ=-10∫ξ=0θφ¯ξ-θdηθϕξdξ,where ηθ=ηθ,0.Letqθ=1,q2,q3,q4Teiω0τ0θ be the eigenvector of A0 belonging to +iω0τ0 and q∗s=V1,q2∗,q3∗,q4∗eiω0τ0s be the eigenvector of A∗0 belonging to -iω0τ0. It is not difficult to show that(34)q2=α21iω0-α22,q3=α32q2iω0-α33-β33e-iτ0ω0,q4=iω0-α11-α12q2-β13e-iτ0ω0q3α14,q2∗=-iω0+α11α21,q3∗=-β13eiτ0ω0+α43q3iω0+α33+β33eiτ0ω0,q4∗=-α14iω0+α44.From (33), we can get(35)q∗s,qθ=V¯1+q2q¯2∗+q3q¯3∗+q4q¯4∗+τ0e-iτ0ω0q3β13+β33q¯3∗.Then we choose(36)V¯=1+q2q¯2∗+q3q¯3∗+q4q¯4∗+τ0e-iτ0ω0q3β13+β33q¯3∗-1.such that q∗,q=1.Next, we compute the coordinates to describe the center manifoldC0 at μ=0. Let xt be the solutions of (31) when μ=0. Define(37)zt=q∗,xt,Wt,θ=xtθ-2Reztqθ.On the center manifold C0, we have(38)Wt,θ=Wz,z¯,θ=W20θz22+W11θzz¯+W02θz¯22+⋯,where z and z¯ are local coordinates for center manifold C0 in the direction of q∗ and q¯∗. We now consider only the real solution xt∈C0 of (31), which gives(39)z˙=iω0τ0z+q¯∗F0,Wz,z¯,0+2Rezqθ=iω0τ0z+gz,z¯,where(40)gz,z¯=q¯∗0F0z,z¯=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋯.Thus,(41)gz,z¯=V¯1,q¯2∗,q¯3∗,q¯4∗F10,xt,F20,xt,0,0T,where(42)F10,xt=-βτ0ϕ10ϕ20,F20,xt=βτ0ϕ10ϕ20.Since(43)xt=Wz,z¯,θ+zqθ+z¯q¯θ,qθ=1,q2,q3,q4Teiω0τ0θ,we have(44)xt=W1t+θW2t+θ00+z1q2q3q4eiω0τ0θ+z¯1q¯2q¯3q¯4e-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+⋯,ϕ40=zq4+z¯q¯4+W2040z22+W1140zz¯+W024z¯22+⋯.From (40)–(43), we have(45)gz,z¯=V¯1,q¯2∗,q¯3∗,q¯4∗×M11z2+M12zz¯+M13z¯2+M14z2z¯M21z2+M22zz¯+M23z¯2+M24z2z¯00+⋯,with(46)M11=-βq2,M12=-βq2+q¯2,M13=-βq¯2,M14=-βW1110q2+12W2010q¯2+W1120+12W2020,M21=βq2,M22=βq2+q¯2,M23=βq¯2,M24=βW1110q2+12W2010q¯2+W1120+12W2020.Therefore,(47)gz,z¯=V¯M11+q¯2∗M21z2+M12+q¯2∗M22zz¯+M13+q¯2∗M23z¯2+M14+q¯2∗M24z2z¯+⋯.Thus, from (40) and (47), we can obtain(48)g20=2τ0V¯βq2q¯2∗-1,g11=τ0V¯q¯2∗-1βq2+q¯2,g02=2τ0V¯βq¯2q¯2∗-1,g21=2βτ0V¯q¯2∗-1βW1110q2+12W2010q¯2+W1120+12W2020.In order to computeg21, we need to compute W20 and W11. From (37)-(40), we have(49)W˙=AW-2Req¯∗0F0qθ,-1≤θ<0,AW-2Req¯∗0F0qθ+F0,θ=0,=AW+Hz,z¯,θ,where(50)Hz,z¯,θ=H20θz22+H11θzz¯2+H02θz¯22+⋯.From (38), (39), (49), and (50), one can obtain(51)2iω0τ0-AW20θ=H20θ,(52)AW11θ=-H11θ.Now for θ∈-1,0,(53)Hz,z¯,θ=-2Req¯∗0F0qθ=-g20qθ+g¯11q¯θz22-g11qθ+g¯02q¯θzz¯+⋯,which on comparing the coefficients with (50) gives(54)H20θ=-g20qθ-g¯02q¯θ,(55)H11θ=-g11qθ-g¯11q¯θ.From (51), (54), and the definition of A, we have(56)W˙20θ=2iω0τ0W20θ+g20qθ+g¯02q¯θ,W˙11θ=g11qθ+g¯11q¯θ.Thus,(57)W20θ=ig20q0τ0ω0eiτ0ω0θ+ig¯02q¯03τ0ω0e-iτ0ω0θ+E1e2iτ0ω0θ,(58)W11θ=-ig11q0τ0ω0eiτ0ω0θ+ig¯11q¯0τ0ω0e-iτ0ω0θ+E2,where E1 and E2 are constant vectors, to be determined. It follows from the definition of A, (54), and (55) that(59)∫-10dηθW20θ=2iω0τ0W200-H200,(60)∫-10dηθW11θ=-H110,where ηθ=η0,θ. From (51) and (52), we obtain(61)H200=-g20qθ-g¯02q¯0+M11,M21,0,0T,(62)H110=-g11qθ-g¯11q¯0+M12,M22,0,0T.Noticing that(63)iω0τ0I-∫-10dηθeiω0τ0q0=0,-iω0τ0I-∫-10dηθe-iω0τ0q¯0=0,and substituting (57) and (61) into (59), we obtain(64)2iω0τ0I-∫-10dηθe2iω0τ0θE1=-βq2βq200.Therefore, we can obtain(65)E1=2iω0-α11-α12-β13e-2iτ0ω0-α14-α212iω0-α22000-α322iω0-α33-β33e-2iτ0ω000-α42-α432iω0-α44-1×-βq2βq200.Similarly, we have(66)E2=-α11α12β13α14α21α22000α32α33+β3300α42α43α44-1×-2βReq22βReq200.Then, we have(67)C10=i2τ0ω0g11g20-2g112-g0223+g212,μ2=-ReC10Reλ′τ0,ρ2=2ReC10,T2=-Im⁡C10+μ2Im⁡λ′τ0τ0ω0.Thus, the properties of the Hopf bifurcation of system (2) can be stated as follows.Theorem 3.
Direction of the Hopf bifurcation is determined byμ2: if μ2>0(μ2<0), then the Hopf bifurcation is supercritical (subcritical). Stability of the bifurcating periodic solutions is determined by ρ2: if ρ2>0(ρ2<0), then the bifurcating periodic solutions are stable (unstable). Period of the bifurcating periodic solutions is determined by T2: if T2>0(T2<0), then the period of the bifurcating periodic solutions increases (decreases).


## 4. Numerical Simulations

In this section, numerical simulations are presented by taking partial parameters from numerical simulations in [14] and they are as follows: A=1, β=0.35, ϖ=0.8, ε=0.15, μ=0.35, γ=0.05, α=0.35, and δ=0.15. Then, the following system can be obtained:(68)dStdt=1-0.35StLt+0.8It-τ-0.15Rt-0.35St,dLtdt=0.35StLt-0.55Lt,dItdt=0.05Lt-0.7It-0.8It-τ,dRtdt=0.15Lt+0.35It-0.2Rt.By virtue of the chosen values of parameters, we can getR0=1.8182 and the unique endemic equilibrium P∗1.5714,0.6712,0.0224,0.5427. Further, we can verify that the conditions indicated in Theorem 2 are satisfied. In this case we can obtain ω0=0.5209, τ0=5.3442, and f′ω02=1.9227≠0.Therefore,P∗1.5714,0.6712,0.0224,0.5427 is asymptotically stable when 0<τ<τ0=5.3442 according to Theorem 2, which can be shown as the numerical simulation in Figure 1. In this case, propagation of the computer viruses can be controlled easily. However, system (68) undergoes a Hopf bifurcation at P∗1.5714,0.6712,0.0224,0.5427 when τ passes through the critical value τ0. This property can be illustrated by the numerical simulation in Figure 2 and propagation of the computer viruses will be out of control. In addition, according to (67), we have μ2=0.8048>0, ρ2=-1.2498<0, and T2=-0.7991<0. Thus, we know that the direction of the Hopf bifurcation at τ0 is supercritical; the period of the bifurcating periodic solutions decreases; and the bifurcating periodic solutions are stable. From the viewpoint of biology, if the periodic solutions bifurcating from the Hopf bifurcation are stable, then the susceptible, exposed, infectious, and recovered nodes in system (68) may coexist in an oscillatory mode. This phenomenon is not welcome in a real network. We regret to say that only supercritical Hopf bifurcation is identified in our numerical case study. However, it should be pointed out that subcritical Hopf bifurcations are quite usual in dynamical systems with time delay, and the existence of an unstable periodic solution makes the dynamics even more intricate, which has been discussed in population dynamics in the literature [22].Figure 1
P ∗ is locally asymptotically stable when τ=5.2482<τ0=5.3442 with initial functions “1.55, 0.79, 0.136, 0.75.”Figure 2
P ∗ lose its stability when τ=6.1752>τ0=5.3442 with initial functions “1.55, 0.79, 0.136, 0.75.”It should be also pointed out that we know that onset of the Hopf bifurcation can be postponed by properly increasing values of the parametersδ, α, and μ based on numerical simulations. In reality, cure rate of the exposed node δ and cure rate of the infectious nodes can be properly enhanced by updating of antivirus software on nodes; disconnecting rate of every node μ can be properly increased by the managers of a real network. Therefore, propagation of the computer viruses in system (2) can be controlled by timely updating of antivirus software on nodes and timely disconnecting nodes from a real network when the connections are unnecessary.


## 5. Conclusions

A delayed computer virus model is proposed based on the literature [14] considering that the outbreak of computer virus usually lags. By theoretical analysis, the critical value of the delay τ0 where the Hopf bifurcation occurs is obtained and our results show that the time delay plays an important role in the stability of the considered computer virus model in the present paper.When the value of the delayτ<τ0, the endemic equilibrium of the model is asymptotically stable and the propagation of the computer virus is divinable. However, when the value of the delay τ>τ0, the endemic equilibrium of the model will lose its stability and a Hopf bifurcation occurs. In this case, propagation of the computer viruses will be out of control. Therefore, we should postpone occurrence of the Hopf bifurcation by taking some effective measures, such as timely updating of antivirus software on nodes and timely disconnecting nodes from the network when the connections are unnecessary. In addition, properties of the Hopf bifurcation are also investigated by using the normal form theory and center manifold theorem.Of course, various factors in a real network can affect propagation of the computer viruses. Among them, time delay is important and the present paper thus concentrates on analyzing it. Other network impacts of computer viruses such as network topology, network eigenvalue, and patch forwarding cannot be neglected, and they will be a major emphasis of our future research. Also, stability of the endemic equilibrium and existence of the periodic solutions are investigated only at the critical valueτ0. As stated in the literature [22], it is an open question yet whether there exists stable endemic equilibrium for even larger time delay. This is the other interesting future research direction for us.



---

*Source: 9360430-2017-06-20.xml*