Impact of Discontinuous Treatments on Disease Dynamics in an Sir Epidemic Model

We consider an SIR epidemic model with discontinuous treatment strategies. Under some reasonable assumptions on the discontinuous treatment function, we are able to determine the basic reproduction number R 0 , confirm the well-posedness of the model, describe the structure of possible equilibria as well as establish the stability/instability of the equilibria. Most interestingly, we find that in the case that an equilibrium is asymptotically stable, the convergence to the equilibrium can actually be achieved in finite time, and we can estimate this time in terms of the model parameters, initial sub-populations and the initial treatment strength. This suggests that from the view point of eliminating the disease from the host population, discontinuous treatment strategies would be superior to continuous ones. The methods we use to obtain the mathematical results are the generalized Lyapunov theory for discontinuous differential equations and some results on non-smooth analysis. 1. Introduction. Infectious diseases remain to be one of the main sources of deaths for the human beings. Advances of knowledge on the existing infectious diseases, including knowledge about the mechanisms of transmissions, can help reduce the spread and hence the deaths. One important approach to understand disease transmission mechanisms in the population level is mathematical modeling, and differential equations play a crucial role in this regards because such equations describe how the rate of changes of sub-populations in host population depends on the main parameters as well as on the sub-populations themselves. There have been many differential equations models and the books [4] and [11] provide good coverages of some basic and classic models.

1. Introduction.Infectious diseases remain to be one of the main sources of deaths for the human beings.Advances of knowledge on the existing infectious diseases, including knowledge about the mechanisms of transmissions, can help reduce the spread and hence the deaths.One important approach to understand disease transmission mechanisms in the population level is mathematical modeling, and differential equations play a crucial role in this regards because such equations describe how the rate of changes of sub-populations in host population depends on the main parameters as well as on the sub-populations themselves.There have been many differential equations models and the books [4] and [11] provide good coverages of some basic and classic models.
When an infectious disease emerges in a host population, it is desirable to consider some control measures.A natural question arises: how would such a measure 98 ZHENYUAN GUO, LIHONG HUANG AND XINGFU ZOU affect the disease dynamics?A natural way to answer this question is to incorporate the control measure into the existing mathematical model that has been well understood and investigate how the measure will affect the behavior of the solutions to the model.Note that a control measure here is in a general sense and can include vaccination, curing treatment, quarantining, and isolation etc.For effects of vaccination on disease dynamics in some models, see, e.g., [1,2,5,9] and the references therein; for some models studying the impact of quarantine and/or isolation, see, e.g., [6,15,16,19,20,22] and their references.For general treatment, Brauer [10] considered a model for a heterogeneous population with a treatment, and investigated the impact of mixing and treatment.Wang [21] added a treatment function to a classic SIR model and observed the backward bifurcation which would not exist in the absence of the treatment.Motivated by [21], Zhang [23] extended the study to an SIS model and also obtained backward bifurcation.This work is motivated by [21] where the treatment function is assumed to be continuous.In reality, the treatment strategy usually is not smooth and even not continuous because, due to limited resources, usually there are some restrictions on treatment strengths.In this paper, we use the model in [21], but we adopt discontinuous treatment function.Since the resulting model is a discontinuous system, we need to make use of the theory for discontinuous differential equations and some results on non-smooth analysis from [7,8,13,14,17].Under some reasonable assumptions on the discontinuous treatment function, we explore the well-posedness (Section 2), and structure as well as stability of equilibria (Section 3).The basic reproduction number to such a discontinuous epidemic model is also calculated and its threshold role is confirmed too.As the most novel part of this work, we find in Section 4 that the convergence to an equilibrium can be achieved in finite time.This is especially significant when the equilibrium is the disease free equilibrium since this provides possibility to drive a disease to extinction in finite time, a conclusion in strong contrast to the results for models with continuous treatment functions where convergence to an equilibrium is in the sense of "asymptotic" (i.e., as t → ∞).Our results even give explicit formula for such a finite time, and this allows us to investigate how the model parameters, initial sub-populations and the initial treatment strength affect this time.This suggests that from the view point of eliminating the disease from the host population, discontinuous treatment strategies would be superior to continuous ones.We point out that the phenomenon of convergence to an equilibrium in finite time has been reported for discontinuous neural network models in [18] which stimulated our work on this topic.

Model description and preliminaries.
Consider the following SIR model with treatment proposed in [21] : where state variable S, I and R denote the sub-populations of susceptible, infective and recovered classes respectively in a host population.Among the model parameters which are all positive constants, A is the recruitment rate of the population; d is the natural death rate of the population; γ is the natural recovery rate of infective individuals; is the disease-related death rate; λ is the infection coefficient.Function h(I) represents the treatment rate, which was assumed to be continuous in I in [21].But here, for realistic considerations as mentioned in the introduction, we allow some possible jump discontinuities by assuming the following (A1) h(I) = ϕ(I)I, where ϕ : R + → R + is non-decreasing and has at most a finite number of jump discontinuities in every compact interval.
Remark 1.Without loss of generality, we can always assume that ϕ is continuous at 0, since otherwise, we can modify the value of ϕ at 0 to be ϕ(0 + ) and this will not change the fact of h(0) = 0, and will bring no change to (1) and ( 2).
Since the first two equations in model ( 1) are independent of the variable R, it is sufficient to consider the following sub-system: Due to the discontinuity on the right hand side in (2), many results in the classical theory of ordinary differential equations can not be applied here.To proceed, firstly, we need to give an appropriate definition of solution for model (2).Here, we adopt the definition of solution in the sense of Filippov [17].
Note that under (A1), co[h(I)] is an interval with non-empty interior when h(I) is discontinuous at I, while co[h(I)] = h(I) is a singleton when h(I) is continuous at I. It is easy to see that the map (S, I) → (A−dS −λSI, λSI −(d+γ + )I −co[h(I)]) is an upper semi-continuous set-valued map with non-empty compact convex values ( [7], p.67, Lemma 1).By the Measurable Selection Theorem ( [7], p.90,Thm.1), we know that if (S, I) is a solution of model ( 2) on [0, T ), where T ∈ (0, +∞], then, there exists a measurable function m(t) ∈ co[h(I(t))] such that Remark 2. It is obvious that a) the measurable function m(t) in ( 4) is uniquely determined by (S(t), I(t)) up to a set of measure zero in [0, T ); and b) m(t) is continuous for all t ∈ [0, T ) if and only if (S(t), I(t)) is continuously differentiable for all t ∈ [0, T ).
For biological reason, we need to prove the positiveness and boundedness of solutions of model (2) with positive initial values.
Proof.By the definition of solution of model (2) in Filippov sense, (S(t), I(t)) is a solution of differential inclusion (3).From the S equation in (3), we have This together with S(0) = S 0 ≥ 0 shows that S(t) ≥ 0, t ∈ [0, T ).Note that (A1) implies that co[h(0)] = {0} and h(I) is continuous at 0. By the continuity of ϕ at I = 0 (see Remark 1), there exists a positive constant δ such that when |I| < δ, ϕ(I) is continuous and the differential inclusion in (3) becomes the following differential equation with continuous right hand side: Now if I 0 = 0, it follows from (5) that I(t) = 0 for all t ∈ [0, T ).If I 0 > 0, we claim that I(t) > 0 for all t ∈ [0, T ).Otherwise, let t 1 = inf{t : I(t) = 0}.Then, t 1 > 0 and I(t 1 ) = 0.It follows from the continuity of I(t) on [0, T ) that there exists a positive constant θ such that t 1 − θ > 0 and 0 which is a contradiction.Therefore, I(t) > 0 for all t ∈ [0, T ).The next result addresses the global existence and boundedness of solutions to the model (2).Proposition 2. Suppose that (A1) is satisfied.Then, for any S 0 ≥ 0 and I 0 ≥ 0 there is at least one solution (S(t), I(t)) to the model (2) satisfying S(0) = S 0 and I(0) = I 0 .Furthermore, any such solution is bounded and exists for all t ∈ [0, +∞).
3. Equilibria and their stability.By an equilibrium of model ( 2), we mean a constant solution of model ( 2), (S(t), I(t)) = (S * , I * ), t ∈ [0, +∞).Clearly, (S * , I * ) is an equilibrium of model ( 2) if and only if Thus, if (S * , I * ) is an equilibrium of model ( 2), then there exists a constant ξ * ∈ co[h(I * )] such that Obviously, such a constant ξ * is unique, being given by ξ Suppose that (A1) holds.In order to obtain the equilibria of model ( 2), we need to solve the following inclusion: Obviously, the disease free equilibrium E 0 = (A/d, 0) always exists.An endemic equilibrium satisfies Solving the first equation of ( 9) for S in terms of I gives S = A/(d + λI).Substituting this into the second equation (inclusion), we have Denote and let which is the basic reproduction number of the model (2).We will see in the sequel that the existence of endemic equilibrium is closely related to the size of R 0 .
We next show that Ĩ is the unique positive solution of (10).Set I * 1 = Ĩ and assume that I * 2 = I * 1 is another positive solution of (10).Then, there exist From the monotonicity of ϕ (see (A1)), it follows that Subtraction of the two equations in (12) results in > 0, which is a contradiction.Therefore, (10) has the unique positive solution Ĩ, and the proof of the lemma is completed.
A direct consequence of Lemma 1 is the following uniqueness theorem for endemic equilibrium.By (A1) and Remark 1, we can analyze the stability of the model at the disease free equilibrium E 0 by investigating the eigenvalues of the Jacobian matrix of ( 2) at E 0 .This matrix is calculated as It is clear that the stability of E 0 is fully determined by the sign of the term λA/d−(d+γ + )−ϕ(0): E 0 is asymptotically stable if λA/d−(d+γ + )−ϕ(0) < 0; it is unstable if λA/d − (d + γ + ) − ϕ(0) > 0. The above stability criteria can be stated in terms of R 0 .
Theorem 2. Assume that (A1) holds.Then the disease free equilibrium E 0 is asymptotically stable if R 0 < 1 and it becomes unstable when R 0 > 1.Now we turn to the stability of the unique endemic equilibrium E * = (S * , I * ) which exists if R 0 > 1.We can show that R 0 > 1 is actually also a necessary condition for the existence of the endemic equilibrium E * = (S * , I * ).Indeed, assuming that E * exists, it follows from (9) that λA = (η * + d + γ + )(d + λI * ), where η * ∈ co[ϕ(I * )].Thus, by the monotonicity of ϕ, we have which implies that R 0 > 1. Therefore R 0 > 1 is a sufficient and necessary condition for the existence of the unique endemic equilibrium E * = (S * , I * ).
Assume that R 0 > 1 and that ϕ is differentiable at I * .Then the Jacobian matrix of (2) at the endemic equilibrium E * can be calculated as Theorem 3. Suppose that (A1) holds.If R 0 > 1 and that ϕ is differentiable at I * , then the endemic equilibrium E * is asymptotically stable.
The above stability results are local.Moreover, for E * it is assumed that ϕ is differentiable at I * .In the sequel, we will show that E 0 is indeed globally asymptotically stable when R 0 ≤ 1; and E * is globally asymptotically stable when R 0 > 1, regardless of whether or not ϕ is differentiable at I * .To this end, we need to apply the Lyapunov theory for discontinuous systems (see, e.g., [7,8]).We begin by introducing a LaSalle-type invariance principle.
Consider a system described by the following differential inclusion where F is an upper semi-continuous set-valued map from R n to R n with compact and convex values.We also assume 0 ∈ F (0), that is, 0 is an equilibrium of (13).
A Lyapunov function for ( 13) is a smooth function V : R n → R satisfying the following conditions: (L1) Positive Definiteness: V (x) > 0 for all x = 0; in addition, V (0) = 0. (L2) Properness: the sublevel set {x ∈ R n : V (x) ≤ a} is bounded for every a ≥ 0; (L3) Strong Infinitesimal Decrease: A set W is said to be weakly invariant for (13) if for any x 0 ∈ W , there is at least one solution x(t) satisfying x(0) = x 0 such that x(t) ∈ W for all t at which x(t) exists.
Let V be a Lyapunov function for (13).For any l > 0, by (L1) and the continuity of V , the level set {x ∈ R n : V (x) ≤ l} contains a neighborhood of 0. Denote by V l the largest connected component of the level set {x ∈ R n : V (x) ≤ l} that contains 0. The following LaSalle-type invariance principle is from Theorem 3 in [8].
Lemma 2. Assume that V : R n → R is a Lyapunov function for (13) and let Denote by M the largest weakly invariant subset of Z V L l .Let x 0 ∈ L l and x(t) be any solution with x(0) = x 0 .Then dist(x(t), M ) → 0 as t → +∞.In particular, if M = {0} and l = +∞, then 0 is globally asymptotically stable for (13).Now, we are in the position to state and prove the following global stability result for the disease free equilibrium.Theorem 4. Suppose that (A1) is satisfied.If R 0 ≤ 1, then the disease free equilibrium E 0 is globally asymptotically stable.
Proof.In order to apply Lemma 2, we shift the disease free equilibrium E 0 to the origin by letting x = S − A d .Then, (3) is transformed to the following form: Let which is obviously a smooth function.It is easy to verify that (L1) and (L2) are satisfied for V 1 .Denote the right hand side of ( 14) by G(x, I), that is, .
The following theorem deals with the global asymptotic stability of the endemic equilibrium E * when R 0 > 1. where Consider the function This is a smooth function with respect to (x, y).It is easy to verify that (L1) and (L2) are satisfied.Denote It is easy to see that the map H(x, y) is an upper semi-continuous set-valued map with nonempty compact convex values.For any v = (v 1 , v 2 ) ∈ H(x, y), there exists a corresponding function The gradient of V 2 (x, y) is given by Thus, The monotonicity of ϕ implies that (η(t) − η * )y ≥ 0, and hence, That is ∇V 2 (x, I) • v ≤ 0 verifying (L3).Thus, V 2 is a Lyapunov function for (15).Letting ∇V 2 (x, y) • v = 0, we can obtain Z V2 = {(0, 0)} {(0, y) : η(t) = η * , y = 0}.If x ≡ 0, then by the first equation of ( 15) we obtain y = 0. Therefore, for any l > 0 the largest weakly invariant subset of Z V2 L l for ( 15) is M = {(0, 0)}.By Lemma 2, (0, 0) is globally asymptotically stable for (15); that is, E * is globally asymptotically stable for (2).The proof is completed.

4.
Global convergence in finite time.One important feature of discontinuous ODE systems that a smooth ODE system can not have, is that convergence to equilibrium in finite time is possible under some conditions.This topic has been particularly explored recently by Forti et al [18] for discontinuous neural network models.Motivated by this, in the sequel, we investigate the possibility of convergence to equilibrium in finite time for our model (2).To this end, we need to apply non-smooth Lyapunov functions as was done in Forti et.al. [18], which requires a generalization of the notion of gradient.
A function f : R n → R is said to be regular at x if the following hold: (i) it is locally Lipschitz near x; (ii) for all direction v ∈ R n , there exists the usual one-sided directional derivative , where is the generalized directional derivative of f at x in the direction of v.A function f is said to be regular in R n , if it is regular at every x ∈ R n .
Let f : R n → R be locally Lipschitz in R n .Then, by Rademacher Theorem, f is differentiable at almost all (a.a.) x ∈ R n in the sense of Lebesgue measure (see, e.g., [14,18]).For such a function, the Clarke Generalized Gradient (see [14]), denoted by ∂f , is defined by where Ω f is the set of measure zero in which the gradient of f is not differentiable.
A more meaningful and desirable situation is the global convergence to the disease free equilibrium E 0 in finite time.Since (A1) assumes continuity of the treatment function h(I) at I = 0, finite time convergence to E 0 is impossible under (A1).Thus, discontinuity is required for h(I) at I = 0, as is stated in the following assumption: (A3) h : R + → R + is non-decreasing and has at most a finite number of jump discontinuities in every compact interval.Moreover, h(0) = 0 and h(I) is discontinuous at I = 0.A typical treatment function satisfying (A3) is the following This corresponds to an immediate response to the occurrence of a disease with the constant r being an effort strength.Under (A3), by (7), we know that (A/d, 0) is the disease free equilibrium of model From ( 4), there exists a measurable function η ∈ co[h(I)] corresponding to (x(t), where Proof.Let V 1 (x, I) be the same Lyapunov function as in the proof of Theorem 4.Then, evaluating V1 (x(t), I(t)) along the system (19) gives By (A3), we know that η(t) ≥ h(0 + ).Note that d + γ + − λA d ≥ 0.Then, we have dV Integrating both sides of the above inequality from 0 to t, we obtain 0 ≤ V 1 (x(t), I(t)) ≤ V 1 (x(0), I(0)) − Ah(0 + ) d t = Q(S 0 , I 0 ) − Ah(0 + ) d t.

Conclusion and discussion.
We have revisited the SIR model with treatment considered by [21].But unlike in [21] where treatment function is assumed to be continuous, our main concern here is the impact of the adoption of a discontinuous treatment function.
Our results on the disease free equilibrium E 0 show that when the basic reproduction number R 0 of the model is less than one, as is expected, the disease free equilibrium is globally asymptotically stable.Note that under assumption (A1), R 0 depends on ϕ(0) by (11), which is decreasing function of ϕ(0) (the initial treatment rate).Thus, a larger initial treatment rate will help eliminate the disease, and the formula (11) determines how large ϕ(0) should be.
Under (A1) and when R 0 > 1, the disease free equilibrium becomes unstable.However, the existence, uniqueness and global asymptotic stability of an endemic equilibrium all follow, regardless of whether I * is a continuous or discontinuous point of ϕ(I).
What we believe is most interesting and most novel in this paper are the results on the convergence to an equilibrium in finite time.This is impossible if a smooth treatment function is adopted.Therefore it presents a true advantage of discontinuous treatments.We are also able to establish an estimation of the precise time (finite) it takes for a solution to settle to the equilibrium.This is particularly important and useful for designing treatment strategies aiming to eliminate the disease in finite time.From the expressions ( 16) and (21), one may easily see how the model parameters as well as the initial values and the initial treatment strength will affect the (finite) time it takes to eradicate the disease.Taking (21) as an example, we find that t * is increasing in the magnitude of the initial infectious population I 0 and decreasing in the initial treatment strength h(0 + ), which are all reasonable and natural.Since most, if not all, existing works on disease models with treatments assume continuous treatment functions, our results here on the finite time convergence suggest that it should be worthwhile to reconsider those models by incorporating discontinuous treatment functions.