Skip to main content
Springer Nature - PMC COVID-19 Collection logoLink to Springer Nature - PMC COVID-19 Collection
. 2020 Jul 13;65(1-2):295–319. doi: 10.1007/s12190-020-01392-x

A novel control set-valued approach with application to epidemic models

Lahoucine Boujallal 1, Mohamed Elhia 2,, Omar Balatif 3
PMCID: PMC7355539  PMID: 32837465

Abstract

This work is considered in the framework of studies dedicated to the control problems, especially in epidemiology where the scientist are concerned to develop effective control strategies to minimize the number of infected individuals. In this paper, we set this problem as an asymptotic target control problem under mixed state-control constraints, for a general class of ordinary differential equations that model the temporal evolution of disease spread. The set of initial data, from which the number of infected people decrease to zero, is generated by a new type of Lyapunov functions defined in the sense of viability theory. The associated controls are provided via selections of adequately designed feedback map. The existence of such selections is improved by using Micheal selection theorem. Finally, an application to the SIRS epidemic model, with numerical simulations, is given to show the efficiency of our approach. To the best of our knowledge, our work is the first one that used a set-valued approach based on the viability theory to deal with an epidemic control problem.

Keywords: Non-linear control systems, Epidemic models, Lyapunov functions, Viability theory, Contingent cone, Selections

Introduction and problem statement

Infectious diseases have marked the history of human societies. Throughout the centuries and the world, they have always been the leading cause of death. One thinks of the great plagues of the Middle Ages which emptied cities of their populations and caused urban civilization to regress (a third of the European population disappeared between the end of the fourteenth and the beginning of the fifteenth century [1]). The viruses imported into America by the Spaniards which decimated the local populations more surely than the fights. The Spanish flu which caused more deaths at the beginning of the twentieth century than the First World War. The Acquired Immunodeficiency Syndrome, better known under its acronym AIDS, which is still the most deadly infectious disease in the world, and which has claimed nearly 25 million victims since its appearance in 1981. We also think of the Coronavirus disease (Covid-19), which is now considered as a global pandemic.

To counter the ravages that infectious diseases can cause, public health decision-makers must have relevant tools to assist them in their decision-making. In this context, mathematicians, epidemiologists and immunologists have been working together for quite some time to create mathematical models that allow competent authorities to prepare in advance to respond quickly and effectively if an outbreak occurs.

The first model was developed by Daniel Bernoulli in 1760 for smallpox [2]. But, the foundations of the epidemic modelling based on the compartmental models were established by Sir Ronald Ross, W. H. Hamer, W. O. Kermack and others. In these models, the population is divided into classes that contain individuals with the same epidemiological status. An interesting overview of the history of mathematical models in epidemiology can be found in [36]. Currently, epidemiological models are increasingly used and present a powerful tool for studying complex systems. Its contribution to the fight against epidemics is indisputable. It enables us to understand the mechanisms of transmission, to study the characteristics of an epidemic, to predict its evolution and to evaluate different intervention strategies to find the best control program.

In the literature, several works have addressed the problems of controlling the spread of infectious diseases based on epidemiological models with compartments. Di Giamberardino and Iacoviello [7] proposed an epidemic control problem based on a three compartments model and a vaccination strategy with a cost index that weights differently the control depending on the severity of the disease. Buonomo et al. [8] focus on an epidemic model which incorporates a non-linear force of infection and two controls: an imperfect preventive vaccine given to susceptible individuals and treatment given to infectious. Work done by Zakary et al. [9], devise a discrete time compartment model depicting the spread of infectious diseases in various geographical regions that are connected by any kind of anthropological movement. The authors used two control variables which represent the effectiveness rates of vaccination and travel-blocking operation, and they focus to control the outbreaks of an epidemic that affects a hypothetical population belonging to a specific region. Another work that concerns the control of the Ebola virus disease was proposed by Mhlanga [10]. The mathematical model includes control functions representing educational campaigns in their respective patches, with one patch having more effective controls than the other. Moualeu et al. [11] presented a deterministic model for the transmission dynamics of Tuberculosis in the context of weak diagnosis capacity. Optimal control theory is used to obtain a cost-effective balance of two different intervention methods. Degang et al. [12] presented a novel SIVRS mathematical model for infectious diseases, which takes into account the virus variation factors. Authors proposed an optimal control problem to maximize the recovered agents with the limited resource allocation. Existence of a solution to the optimal control problem is given based on Pontryagin’s Minimum Principle. Bolzoni et al. [13] proposed a problem of minimizing epidemic size and duration in SIR models. Optimal control through either vaccination or isolation is investigated by application of Pontryagin’s Maximum Principle. Using a basic SEI model with saturated incidence rate, Baba et al. [14] studied the effect of optimal controller and awareness on the dynamic of Tuberculosis. In these works and in many others (see [1520] and the references therein) the authors have solved the formulated control problems using classical results of control theory, in particular the famous Pontryagin’s minimum principle [21].

In this paper, we propose an alternative approach which is direct and allows characterizing one or more controls aimed at achieving the null number of infected people at a final time. Our approach is based on viability theory [22], which allows the adaptation of the evolution of a dynamic system to the restrictions imposed on the state and the control. We are interested here in compartment models that represent the evolution of an infectious disease in a population, which are written as follows

x˙=f(x,y,u), 1a
y˙=ψ(x,y,u), 1b

with initial conditions

xi(0)=xi0for eachi=1,,nandy(0)=y0, 1c

for an integer n and where the xi’s denotes the number of individuals in a population such as susceptible and recovered, while y denotes the number of infected individuals. Therefore, the state (xy) evolves in the subset

ΩR+n×R+,

and the control u takes values within constraints subset

Ui=1p[0,uimax], 2

where uimax are positive numbers and p an integer. Both functions f and ψ map Rn×R×Rp into Rn and R respectively, and they are sufficiently smooth. Then, our control problem is stated as follows

Problem: For all (x0,y0)ΓΩ, find a control u¯ such that:

u¯:[0,)U, 3a
(x¯(t),y¯(t))Ω,for allt0, 3b
limty¯(t)=0, 3c

where (x¯,y¯) denotes a solution of system (1) for control u¯.

The problem (3) has been treated by Kassara [23] and Kassara and Moustafid [24] in the framework of the set-valued approach for a class of ODE immunotherapy models. Such models was expressed by a non-linear control system (1), in the particular case, where the dynamics f (1a) is affine in control u, and the Eq.(1b) is given as follows

y˙=yψ(x,y). 4

Although the approach proposed by the authors can be used for several models of cancer treatment, it is not appropriate for a large class of models where the state y can be controlled to zero by acting on its dynamics directly, especially for epidemic models, see for instance [25, 26].

Inspired by the work done in [27], we provide a unified set-valued approach in order to deal with the problem (3) for a more general class of ODE systems, where f can be non-affine with respect to the control and the dynamics ψ depend on the control variable. More precisely, our contribution here is to seek conditions under which, any Lyapunov function φ (defined independently from dynamics ψ, in a sense specified below) can generate a domain Ωφ and a closed-convex valued multifunction Gφ in such a way the controls that solves the problem (3) from Ωφ are given via selections of multifunction Gφ.

Such Lyapunov functions consist of C1 real-valued functions φ defined on R+×R and satisfying

h:[0)R+is differentiable andφ(h(t),h˙(t))0for allt0,limth(t)=0. 5

Note such functions φ depend only upon subset R+, we thereby call them R+-Lyapunov functions. This type of functions was firstly introduced in [27] dealing with constrained asymptotic null-controllability.

Our approach takes advantage of the unified framework of viability theory [28] and set-valued analysis [29, 30]. The main advantage is that for any initial data starting from Ωφ the state stays in this region overcoming the concern mentioned in (3c). This is due to the fact that, the multifunction of regulation Gφ includes the tangential condition for viability. Even when initial data start outside the region Ωφ it is possible to reach it in a certain instant. Moreover, it is of interest to stress that in the case of convex constraints, the multifunction Gφ bring us an explicit continuous selections and universal formula, notably the minimal selection despite its discontinuity. Another advantage is that this set-valued framework can be considered as a unified setting for dealing with problem (3) for an important class of ODE models including, for instance, both cancer and epidemic models.

Throughout the paper, the inner product on Euclidean space is denoted Inline graphic, and corresponding norm Inline graphic. For a vector z we denote by zi its ith component. Let T be a linear operator and we denote its adjoin operator by T and its norm by T.

Also, we consider the notation

xψψx1,,ψxnanduψψu1,,ψup.

The remainder of the paper is structured as follows: In Sect. 2 we provide some notations and preliminary lemmas. Section 3 will be devoted to presenting our control set-valued approach, then Sect. 4 to prove the existence of continuous selections under sufficient conditions. To show the efficiency of our approach we propose, in Sect. 5, an application to a controlled SIRS epidemic model. We conclude this paper in Sect. 6.

Definitions and preliminary results

In this section, we present some viability concepts and mathematical tools to be used next.

Let K be a nonempty subset of an Euclidean space RN. The contingent cone to subset K at point xK is defined by

TK(x)={yRNlim infε0d(x+εy,K)ε=0},

where d(y,K)inf{|z-y||zK}. It is useful to note the following properties of the contingent cones:

  1. If xint(K), the interior of K, then TK(x)=RN.

  2. TK(x) is closed for each xK.

  3. If K is convex, then TK(x) is convex for each xK.

  4. If K is closed and convex then the map TK(·) from K to subsets of RN is lsc.

In (d) above, the short-term lsc means that multifunction TK(·) is lower semi-continuous, i.e. for every xK and any sequence of elements xk of K converging to x, it holds that:

for each yTK(x), there exist a sequence ykTK(xk) that converges to y.

Next we provide a result which is fundamental for building our control strategy. In [29] the authors gave a useful characterization of the contingent cone in case of subsets defined by inequalities:

Lemma 1

Let ϕ:RNR be a differentiable mapping for an integer N. Given a closed convex subset K of RN and let D{zK|ϕ(z)0} and z0D. Suppose that there exists ξ0TK(z0) such that dϕ(z0)ξ0<0, then

ξTD(z0)ξTK(z0)anddϕ(z0)ξ0ifϕ(z0)=0, 6

where dϕ(·) denotes the differential operator of ϕ.

Let ξ:RNRN and consider the following system:

z˙=ξ(z),z(0)=z0. 7

The subset D is said to be locally viable under system (7), if for all z0D there exist t¯>0 and a solution to system (7), z¯(·) on [0,t¯) such that z¯(t)D for all t[0,t¯).

The following lemma gives us sufficient conditions, in terms of contingent cones, for which the property above holds true:

Lemma 2

Assume that function ξ is continuous on the closed subset D. Then D is locally viable under system (7) if and only if the following tangential condition holds,

ξ(z)TD(z)for eachzD. 8

The global solution of the system (7) is provided by using the notion of linear growth. A smooth function ξ mapping a subset D of RN into RN is said to be of linear growth on subset D, if there exists a constant k>0 such that

|ξ(z)|k(1+|z|)for allzD.

However, such a condition is not required when the subset D is bounded, or when the function ξ is bounded on this subset.

By selection of a multifunction Q it is meant a function ξ such that ξ(z)Q(z) for all x. We cite at this opportunity the famous Michael selection theorem [30] which will be repeatedly used in the paper:

Lemma 3

If a multifunction Q is lsc and has closed convex values, then for all (z0,v0) such that v0Q(z0), there exists a continuous selection σ of Q which satisfies σ(z0)=v0.

The minimal selection (well defined if Q has closed convex values) of the map Q is given by

ξmin(z)=πQ(z)(0)for allzD. 9

Here π stands for the operator of best approximation that is defined on K as follows:

|x-πK(x)|=d(x,K)for allxRn.

It is of interest to notice that the minimal selection is rarely continuous, see [30].

Next, given the control system

z˙=ξ(z)+G(z)v,vV(z), 10

where v takes values in Rm and denotes the control, G:RNL(Rm,RN), and V(·) stands for the set-valued map of constraints, define the feedback map

G(z)vV(z)|ξ(z)+G(z)vTD(z). 11

Assume that ξ and G(·) are continuous on D, then any continuous selection of the feedback map G provides a control law that leads to a viable solution to system (10) in D. This is so for minimal selection whenever it is continuous. Otherwise we may use [28, Theorem 4.3.2] as follows:

Lemma 4

Assume that the feedback map G is lsc with closed convex values. Then system (10) with feedback control v=πG(z)(0) has a locally viable solution in D.

A control set-valued approach

In this section, we restate the problem (3) in the context of viability theory and set-valued analysis. Throughout the paper the control system (1) is assumed to satisfy the permanently hypotheses below:

Functionsfandψare continuous and oflinear growth on constraint subsetΩ×U.

For the sake of conciseness, the elements of our control strategy proceeds as follows:

  • (i)

    Let φ be an R+-Lyapunov function in the sense of (5).

  • (ii)
    Set
    Dφ(x,y,u)Ω×U|φ(y,ψ(x,y,u))0, 12
    and
    Ωφπ1(Dφ), 13
    where π1 denotes the mapping (x,y,u)(x,y).
  • (iii)
    Define the feedback map given for each (x,y,u)Dφ as follows
    Gφ(x,y,u)vV|(f(x,y,u),ψ(x,y,u),v-αu)TDφ(x,y,u), 14
    where the control v has values in closed subset
    Vj=1p[0,αujmax],
    with α stands for a positive number.
  • (iv)

    Pick up a continuous selection of the map Gφ.

Consequently, we are ready to prove the following result:

Theorem 1

Any continuous selection of the map Gφ provides a solution of the problem (3) for each (x0,y0)Ωφ.

Proof

Let g be such a selection. Then, for all(x,y,u)Dφ, the following tangential condition holds

(f(x,y,u),ψ(x,y,u),g(x,y,u)-αu)TDφ(x,y,u).

According to lemma (2), the following system

x˙=f(x,y,u),x(0)=x0y˙=ψ(x,y,u),y(0)=y0,u˙=g(x,u)-αu,u(0)=u0, 15

has a solution (x¯,y¯,u¯) which is viable in Dφ, for each (x0,y0,u0)Dφ. This solution is global as both f and ψ possess linear growth and g is bounded. It follows that y¯ has values in R+ and satisfies φ(y¯(t),y¯˙(t))0 for all t. Then, thanks to (5), y¯(t)0 when t, and thereby u¯ and (x¯,y¯) satisfy condition (3a) and (3b) respectively.

In consequence, we are led to design the feedback map Gφ(·) of Eq. (14). For that purpose we need first to compute the contingent cone TDφ(·). In all what follow, we assume that both

dynamicsψand Lyapunovfunctionφare differentiable.

Set then for all (x,y,u)Dφ

ϕ(x,y,u)φ(y,ψ(x,y,u)). 16

It follows that subset Dφ given by (12) can be written as

Dφ(x,y,u)Ω×U|ϕ(x,y,u)0. 17

Using [29, Chapter 4] on tangent cones, we get the contingent cone of subset Ω×U as follows

TΩ×U(x,y,u)=TΩ(x,y)×TU(u), 18

for all (x,y,u)Ω×U, where

(θ,ρ)TΩ(x,y)θi0ifxi=0fori=1,,n,ρ0ify=0, 19

and

ξTU(u)ξi0ifui=0fori=1,,p,ξi0ifui=uimaxfori=1,,p. 20

Moreover, we can easily see that the qualification condition in Lemma 1 is satisfied whenever the following conditions hold:

for all(x,y,u)Dφ,such thatϕ(x,y,u)=0,there exists(θ,ρ,ξ)TΩ×U(x,y,u)such that:xϕ(x,y,u),θ+ρϕy(x,y,u)+uϕ(x,y,u),ξ<0. 21

In the following result, we provide an expression of the contingent cone TDφ(·).

Lemma 5

Suppose that condition (21) is satisfied. Then, for each (x,y,u)Dφ, we have

(θ,ρ,ξ)TDφ(x,y,u)(θ,ρ)TΩ(x,y),ξTU(u),and ifϕ(x,y,u)=0xϕ(x,y,u),θ+ρϕy(x,y,u)+uϕ(x,y,u),ξ0. 22

Proof

Thanks to Lemma 1 and condition (21), we can get easily the expression (22).

Note that the partial differentials of ϕ are given on Dφ by

xϕ(x,y,u)=xψ(x,y,u)φz(y,ψ(x,y,u)), 23
ϕy(x,y,u)=φy(y,ψ(x,y,u))+ψy(x,y,u)φz(y,ψ(x,y,u)), 24

and

uϕ(x,y,u)=uψ(x,y,u)φz(y,ψ(x,y,u)). 25

Subsequently, we need to define the following functions and maps,

φ(x,y,u)φz(y,ψ(x,y,u))xψ(x,y,u),f(x,y,u)+φy(y,ψ(x,y,u))+ψy(x,y,u)φz(y,ψ(x,y,u))×ψ(x,y,u)-φz(y,ψ(x,y,u))αu,uψ(x,y,u), 26
mφ(x,y,u)uψ(x,y,u)φz(y,ψ(x,y,u)), 27

and

Cφ(x,y,u){vV|φ(x,y,u)+mφ(x,y,u),v0}, 28

for all (x,y,u)Ω×U. In addition, we consider the assumption

xi=0andy0fi(x,y,u)0for alluU,y=0andx0ψ(x,y,u)0for alluU 29

for i=1,,n.

Then, we are ready to determine a useful expression of the feedback map Gφ(·) given by (14).

Lemma 6

Under condition (21) and (29), we have for all (x,y,u)Dφ

Gφ(x,y,u)=Vifϕ(x,y,u)<0,Cφ(x,y,u)ifϕ(x,y,u)=0. 30

Proof

According to (14) and Lemma 5, we get

vGφ(x,y,u)(f(x,y,u),ψ(x,y,u)TΩ(x,y),v-αuTU(u)and, ifϕ(x,y,u)=0thenxϕ(x,y,u),f+ψϕy(x,y,u)+uϕ(x,y,u),v-αu0, 31

for each (x,y,u)Dφ. Thanks to (19), (20) and (29), it follows that

vGφ(x,y,u)vVand, ifϕ(x,y,u)=0thenxϕ(x,y,u),f+ψϕy(x,y,u)+uϕ(x,y,u),v-αu0. 32

By considering (26), (27) and (28), we can see easily that the last expression is equivalent to (30).

In what follow, we provide a result for which the problem (3) has a solution from the whole domain Ω. For this end, we need to introduce the following map, for each β>0.

Cφβ(x,y,u)vV|φ(x,y,u)+mφ(x,y,u),v-β 33

for each (x,y,u)Ω×U.

Theorem 2

Assume that, for some β>0, the map Cφβ given by (33) has a continuous selection, then problem (3) has a solution for each (x0,y0)Ω.

Proof

Let g be such a selection of the map Cφβ. As condition (21) is satisfied, this map has values included in CφGφ. Then the selection g is also a continuous selection of Gφ. By Theorem 1, Problem (3) has therefore solution from Ωφ. The rest of the proof is devoted to show that problem (3) has solution from Ω\Ωφ.

Now, let (x0,y0) belong to Ω\Ωφ. Thereby ϕ(x0,y0,u)>0, for all u such that (x0,y0,u)Ω×U. Condition (29) implies that (f(x,y,u),ψ(x,y,u))TΩ(x,y) for all (x,y)Ω. Since g(xyu) belongs to V, then g(x,y,u)-αuTU(u). It follows that system

x˙=f(x,y,u),x(0)=x0,y˙=ψ(x,y,u),y(0)=u0,u˙=g(x,y,u)-αu,u(0)=u0,

admits a Ω×U-viable solution (x¯,y¯,u¯) (on horizon 0, as the couple (f,ψ) possess linear growth and g is bounded), for all u0 such that (x0,y0,u0)Ω×U. Let such u0 be given, we readily have

ddtϕ(x¯(t),y¯(t),u¯(t))=xϕ(x¯(t),y¯(t),u¯(t)),x¯˙(t))+y¯˙(t)ϕy(x¯(t),y¯(t),u¯(t))+u¯˙(t),uϕ(x¯(t),y¯(t),u¯(t)),

which yields, according to Eqs. (26) and (27),

ddtϕ(x¯(t),y¯(t),u¯(t))=φ(x¯(t),y¯(t),u¯(t))+u¯˙(t),mφ(x¯(t),y¯(t),u¯(t)).

Whence

ϕ(x¯(t1),y¯(t1),u¯(t1))=ϕ(x0,y0,u0)+0t1[φ(x¯(t),y¯(t),u¯(t))+g(x¯(t),y¯(t),u¯(t)),mφ(x¯(t),y¯(t),u¯(t))]dt.

Since g is also a continuous selection of Cφβ, then

ϕ(x¯(t1),y¯(t1),u¯(t1))ϕ(x0,y0,u0)-βt1.

Thereby

(x¯(t1),x¯(t1),u¯(t1))Dφfort1ϕ(x0,y0,u0)β.

Let x1x¯(t1), y1y¯(t1) and u1u¯(t1).

As a result, the use of Theorem 1 implies that system

s˙=f(s,q,w),s(t1)=x1,q˙=ψ(s,q,w),q(t1)=y1,w˙=g(s,q,w)-αw,w(t1)=u1,

has a solution (s¯,q¯,w¯) on horizon t1,, which ranges in Ω×U and satisfies

q¯(t)0whent.

Consequently, the control given by u¯ on [0,t1] and w¯ on (t1,), achieves the problem stated in (3) from {(x0,y0)}.

On the existence of continuous selection

This section involves the sufficient conditions we need in order to prove the existence of continuous selection of the map Cφβ as required in Theorem 2. For that purpose we will mainly use the Micheal selection theorem cited in section 2 and we assume that

The dynamicsψisC1and haslinear growth on subsetΩ×U.

Definition 1

Let β0. We say that R+-Lyapunov function φ belongs to subset Λβ whenever, for all (x,y,u)Dφ, there exists vV which satisfies the following statement

φ(x,y,u)+v,mφ(x,y,u)<-β, 34

where functions φ and mφ are respectively given in (26) and (27).

Then we can state the following result.

Lemma 7

Let β0 and assume that φΛβ, then both maps Cφβ and Gφ, as given, respectively, by (33) and (30) are lsc with closed-convex valued. (note Cφ0=Cφ).

Proof

It is obvious that Cφβ and Gφ are closed convex valued.

To prove that Cφβ is lsc, we will rewrite it in the context of [29, Proposition 1.5.2] in the form

Cφβ(x,y,u)={vF¯(x,y,u)|f¯(x,y,u,v)G¯(x,y,u)},

where for each (x,y,u)Dφ, we have

F¯(x,y,u)=V,f¯(x,y,u,v)=mφ(x,y,u),v,andG¯(x,y,u)=-,-φ(x,y,u)-β.

Therefore, we can easily verify the hypotheses of the cited proposition as follows:

  1. The map F¯ is lsc with convex values.

  2. f¯ is continuous.

  3. For all (x,y,u)Dφ, the mapping vf¯(x,y,u,v) is affine.

  4. For all (xyu), G¯(x,y,u) is convex and its interior is nonempty.

  5. The graph of the map (x,y,u)Dφint(G¯(x,y,u)) is open.

  6. For all (x,y,u)Dφ, there exists vF¯(x,y,u) such that f¯(x,u,v)int(G¯(x,y,u)).

Thanks to assumption (34), the condition (6) is satisfied. Then the map Cφβ is lsc.

Now we prove that Gφ is lsc.

Let (xn,yn,un)n be a sequence of Dφ that converge to (x,y,u)Dφ and vGφ(x,y,u). We have to look for a sequence (vn)n that satisfies

vnGφ(xn,yn,un)for eachn,andvnv. 35

Suppose that ϕ(x,y,u)<0. Since the function ϕ is continuous and (xn,yn,un)(x,y,u) we can consider the smallest number n0 such that

ϕ(xn,yn,un)<0for allnn0.

Then the sequence defined by

vn=vifnn0,wnifn<n0,

where

wnCφ(xn,yn,un)for alln<n0,

merely satisfies (35) due to the fact that ϕ(xn,yn,un)=0 whenever n<n0. Now assume that ϕ(x,y,u)=0, then vCφ(x,y,u). Since the map Cφ is lsc, it follows that there exists a sequence (vn)n such that vnCφ(xn,yn,un) for each n and vnv. According to (30) we get vnGφ(xn,yn,un) for all n, as required in (35).

Lemma 8

Let β0 and φ belong to Λβ, then the minimal selection of the map Cφβ is continuous.

Proof

By Lemma 7 the map Cφβ is lsc. Then we can use [31, Theorem 4.1] by verifying that the subset

Sε={(x,y,u)Dφ|vCφβ(x,y,u)s.t.|v|Rpε}

is closed in Dφ for all ε>0. Indeed, let ((xn,yn,un))n be a sequence in Sε which converges to (x¯,y¯,u¯). Then there exists a sequence (vn)nRp such that

|vn|Rpεandφ(xn,yn,un)+mφ(xn,yn,un),vn-β 36

for all n. Now, as (vn)n is bounded it has a subsequence (vm)m which converges to v¯. Then by noting that φ and mφ are continuous and letting m in (36), we get

|v¯|Rpεandφ(x¯,y¯,u¯)+mφ(x¯,y¯,u¯),v¯-β.

This implies that (x¯,y¯,u¯)Sε and therefore Sε is closed.

Then we are in a position to state and prove the following result.

Theorem 3

Let β0 and φ belong to Λβ. If β=0 then problem (3) has solution from Ωφ, or else it has solution from the whole domain Ω.

Proof

Suppose that β=0 and Let φΛ0, then, by virtue of Lemma 7, the feedback map Gφ is lsc on Dφ and has a closed convex valued. Furthermore, Lemma 3 implies that Gφ has a continuous selection g. We are now able to use Theorem 1 to conclude that problem (3) has solution from Ωφ.

Now assume that β>0, then, the map Cφβ is lsc and has a closed convex valued (By Lemma 7). As a result the Michael selection theorem (as stated in Lemma 3) yields a continuous selection g of the map Cφβ. The proof ends by using Theorem 2.

The minimal selection of the map Gφ of (30) is given for all (x,y,u)Dφ by the expression

gφ(x,y,u)πGφ(x,y,u)(0)=0ifϕ(x,y,u)<0,πCφ(x,y,u)(0)ifϕ(x,y,u)=0, 37

where the map Cφ is given by Eq. (28). Although mapping gφ is discontinuous, we will see in the next result that it can lead to a slow control law which solve the problem (3).

Theorem 4

Let φ belong to Λ0, then for all (x0,y0,u0)Dφ, system

x˙=f(x,y,u),x(0)=x0,y˙=ψ(x,y,u),y(0)=u0,u˙=gφ(x,y,u)-αu,u(0)=u0, 38

has a solution (x¯,y¯,u¯):[0,)Ω×U which satisfies y¯(t)0 at infinity.

Proof

By Lemma (7), the map Gφ is lsc with closed-convex valued. Then using Lemma (4), system (38) has a solution (x¯,y¯,u¯) over a bounded horizon. Since the couple (f,ψ) has a linear growth and gφ is bounded, then the solution can be extended to infinite horizon. Moreover, since φ is an R+-Lyapunov function, then y¯ goes to zero at infinity.

Application to a controlled SIRS epidemic model

The aim of this section is to apply the approach developed in the previous sections to an epidemic model that includes control terms. We consider a deterministic Susceptible-Infected-Recovered-(re)Susceptible (SIRS) model, where the total population, denoted N, is divided into three compartments:

  • Susceptible (S): healthy individuals who can become infected as a result of their interactions with infected individuals;

  • Infectious (I): individuals who are infected with the disease and are able of transmitting the infection to susceptible individuals;

  • Recovered (R): individuals who were previously infected and recover.

The SIRS model is widely used to model several diseases, such as influenza and malaria, which confer temporary immunity; the recovered individuals lose immunity after a while and become susceptible again. This model is formulated based on the following assumptions

  • The disease incubation is negligible, in this case, once infected, each susceptible individual becomes infectious instantaneously;

  • All recruitment is into the susceptible compartment and occur at a constant rate Λ;

  • The natural death rate, denoted a, is constant across all compartments;

  • The disease is assumed to be fatal for infectious individuals, that is why we define additional death rate c;

  • The transmission of the disease occurs following an adequate contact between a susceptible and infectious. Due to the non-linear contact dynamics in the population, we use the incidence function bSIN to indicate successful transmission of the disease, where b denotes the effective contact rate with infectious individuals in compartment I;

  • The rate of loss of immunity is e;

  • The population size is constant.

The dynamics of our model is governed by the following non-linear system

S˙=Λ-bSIN-aS+eR,I˙=bSIN-(a+c+d)I,R˙=dI-(a+e)R, 39

with S(0)=S0,I(0)=I0 and R(0)=R0.

It is well known that treatment, education and awareness campaign are an important and effective method to control and prevent the spread of various epidemics. Hence, we investigate the impact of these control measures on the spread of an infectious disease by introducing in the system (39) two controls u1 and u2 such as

  • u1: represents the effort of preventing susceptible individuals from becoming infectious individuals. This effort includes awareness program, isolation and any other distancing measurement that can limit contacts between susceptible and infectious people;

  • u2: represents screening and treatment efforts.

So our model with controls is given by

S˙=Λ-(1-u1)bSIN-aS+eR,I˙=(1-u1)bSIN-(a+c+d)I-u2I,R˙=dI-(a+e)R+u2I, 40

with

S0,I0,R0and0uiuimax1fori=1,2.

Let x=(x1,x2)(S,R), yI and u=(u1,u2). We can see that the model (40) satisfies the general form given by system (1), with the dynamic (1a) given as follows

f(x,y,u)Λ-(1-u1)bx1yN-ax1+ex2dy-(a+e)x2+u2y 41

and the infectious dynamic (1b) is given by

ψ(x,y,u)(1-u1)bx1yN-(a+c+d+u2)y, 42

for each (x,y,u)Ω×U, where

Ω=R+2×R+andU=i=12[0,uimax].

It can be easy to see that both functions f and ψ are differentiable. In addition, the condition (29), which gives us the positivity of the model’s solutions, is also satisfied. On the other hand, the global solution exist since both dynamics f and ψ are bounded. In this case, the condition of linear growth is not required.

A convenient R+-Lyapunov function φ, defined in (5) and which satisfies condition (21), can be given by

φ(y,z)z+μy,for all(y,z)R+×R,with(μ>0).

Indeed, if a C1 function h:R+R+ satisfies: φ(h(t),h˙(t))0 for all t, then

h˙(t)+μh(t)0,for allt0.

Using a simple Gronwall inequality, we get h(t)0 at infinity. Now let us express the function ϕ of (16) we get

ϕ(x,y,u)=ψ(x,y,u)+μy,for all(x,y,u)Ω×U,

and let Dφ(x,y,u)|ϕ(x,y,u)0 and Ωφ=π1(Dφ).

To show the efficiency of our approach we consider three scenarios. For each scenario and in order to compute the regulation map Cφβ as expressed in Eq. (33) we first give expressions of the functions φ and mφ, given respectively by Eq. (26) and (27).

All simulations are performed using MATLAB with parameter values given in Table 1. The numerical results are obtained for several initial conditions corresponding to both cases where (x0,y0)Ω\Ωφ and (x0,y0)Ωφ. It should be noted that (x0,y0)Ω\Ωφ refers to when the disease starts to spread in the population, a large part of the population is healthy and the number of infections has not yet peaked. (x0,y0)Ωφ is devoted to describing an advanced stage of the epidemic where the number of infected people is quite large. For the first case, where (x0,y0)Ω\Ωφ, we use the five initial conditions (3500, 5500, 11, 000), (3000, 5000, 12, 000), (2500, 4500, 13, 000), (2000, 4000, 14, 000) and (1500, 3500, 15, 000). In the second case, the initial values are chosen as follows: (3500, 5500, 11, 000), (3000, 5000, 12, 000), (2500,4500,13,000), (2000, 4000, 14, 000) and (1500, 3500, 15, 000).

Table 1.

Parameters description and values

Parameter Description Value References
a Natural death rate 0.00946 [32]
N Total population size 20, 000 Assumed
Λ Recruitment rate aN Assumed
b Transmission rate 0.19 [32]
c Disease induced death rate 0.0353 [32]
d Recovery rate 0.0447 [32]
e Lose of immunity rate 0.23 [33]

First scenario: prevention

Given the major role of contact in transmitting infectious disease from infectious individuals to susceptible people and the importance of prevention programs in limiting the number of new cases, we propose a control strategy based on the control u1, whereas the control u2 is set to zeros. In this case, the expressions of functions φ and mφ are given by

φ(x,y,u)=(1-u1)byNΛ-ax1-(1-u1)bx1yN+ex2+y(1-u1)bx1N-(a+c+d)(1-u1)bx1N-(a+c+d)+μ)+αu1bx1yN, 43
mφ(x,y,u)=-bx1yN. 44

Then the continuous selection, of the map Cφβ, that provides a solution of the problem (3) can be expressed as

g(x,y,u)=minαu1max,max0,-φ+βmφ 45

with β>0 if (x0,y0)Ω\Ωφ and β=0 otherwise.

In Figs. 1 and 2 we depict the evolution over time of the susceptible individuals S, the infected individuals I and the function ϕ in the uncontrolled case and when the control u1 is considered. We use several initial values and we consider both cases (x0,y0)Ω\Ωφ and (x0,y0)Ωφ. In Fig. 1a, it is observed that the number of S decreases sharply when there is no control, while in presence of the control u1 we notice that there is a decrease only during the first ten days and then the number of the susceptible individuals starts to increase. Figure 1b shows a constant increase in the number of infected persons during the first few days followed by a slight decrease when there is no control, but at the end of the control period we can see a clear difference between the case with and without control. In Fig. 2a, b, it can be seen that the prevention control u1 significantly reduces the number of infected people and increases the number of healthy people. The solutions show a decrease in the number of infections of up to 96% at the end of the control period. Also, it should be noted that the control u1 has a great impact on the dynamics of the function ϕ. Indeed, when (x0,y0)Ω\Ωφ the control u1 leads the function ϕ to be negative from a certain time t1. Otherwise, when (x0,y0)Ωφ, the function ϕ remains negative during the full control period.

Fig. 1.

Fig. 1

Cases when the control u1 is applied alone with five initial conditions and (x0,y0)Ω\Ωφ. a Susceptible individuals (S). b Infected individuals (I). c The function ϕ

Fig. 2.

Fig. 2

Case when the control u1 is applied alone with five initial conditions and (x0,y0)Ωφ. a Susceptible individuals (S). b Infected individuals (I). c The function ϕ

Second scenario: screening and treatment

In this case, we propose to provide a treatment for the infected people and facilitate the access to therapeutic programs. Thus, in this scenario we set u1 to zero and we consider only the control u2. The expressions of functions φ and mφ are given by

φ(x,y,u)=byNΛ-ax1-bx1yN+ex2 46
+ybx1N-(a+c+d)bx1N-(a+c+d)+μ+αu2y,mφ(x,y,u)=-y. 47

Thereby the continuous selection of the map Cφβ is given by

g(x,y,u)=minαu2max,max0,-φ+βmφ 48

with β>0 if (x0,y0)Ω\Ωφ and β=0 otherwise.

Figures 3 and 4 display the solutions of S, I and ϕ, without control and when only the control u2 is applied. Compared to the uncontrolled case, there is a steady decrease in the number of the infected individuals. At the end of the control period, the number of I decreases by 100%. With this scenario, we observe a clear decrease in the number of infections from the beginning of the control program whereas in the first scenario there is certainly a regression in the number of infections but more slowly. So in terms of reducing the number of I, the control u2 is a more effective way of eliminating the disease than u1. Also, it should be noted that the control u2 has the same effect on the function ϕ as the control u1.

Fig. 3.

Fig. 3

Cases when the control u2 is applied alone with five initial conditions and (x0,y0)Ω\Ωφ. a Susceptible individuals (S). b Infected individuals (I). c The function ϕ

Fig. 4.

Fig. 4

Case when the control u2 is applied alone with five initial conditions and (x0,y0)Ωφ. a Susceptible individuals (S). b Infected individuals (I). c The function ϕ

Third scenario: combining prevention and treatment

Let us now investigate the effect of combining all controls. With this approach none of the controls is set to zero. In this case, the expressions of functions φ and mφ are given by

φ(x,y,u)=(1-u1)byNΛ-ax1-(1-u1)bx1yN+ex2 49
+y(1-u1)bx1N-(a+c+d+u2)(1-u1)bx1N-(a+c+d+u2)+μ)+α1u1bx1yN+α2u2y,mφ(x,y,u)=mφ,1mφ,2=-bx1yN-y. 50

Then the continuous selection, of the Cφβ, that provides a solution of the problem (3) can be expressed as

g(x,y,u)=g1(x,y,u)g2(x,y,u) 51

where,

g1(x,y,u)=minα1u1max,max0,-φ+β2mφ,1,g2(x,y,u)=minα2u2max,max0,-φ+β2mφ,2,

with β>0 if (x0,y0)Ω\Ωφ and β=0 otherwise.

With this scenario, our results show that combining treatment and prevention measurements can lead to a decrease in the number of infected individuals, an illustrative example is given in Fig.  5a, b. Here, the major advantage of this combination is that less treatment effort is required when it is accompanied by preventive measures. One can observe a significant difference in u2 when it is used alone (Fig. 6a) and when it is used with u1 (Fig. 6b). It can be inferred that in some cases where screening and treating infectious individuals can be more expensive, good preventive measurements with less effort of treatment can also reduce the number of infected people.

Fig. 5.

Fig. 5

Infected individuals without control and when controls u1 and u2 are all used. a Case (x0,y0)Ω\Ωφ. b Case (x0,y0)Ωφ

Fig. 6.

Fig. 6

The control u2 when used alone and when it is coupled with u1. a The control u2 when used alone. b The control u2 when combined with u1

Finally, the results obtained in this application to the SIRS model show the effectiveness of our approach. In the three proposed scenarios, we note that the control terms obtained, allow the number of infected people to tend to zero at a finite time.

Conclusion and future research

This work contributes to a growing literature on controlling epidemics spread. We have proposed a unified approach based on viability theory and set-valued analysis to deal with a control problem of a general class of epidemiological models. We have established theoretical results that show the existence and give the expression of continuous selections used to characterize our controls. The great interest of our approach is the simplicity of deriving the explicit formulas of the controls via continuous selections of adequately designed regulation map, unlike to the optimal control approach, which requires showing the existence of the optimal solution (control and state), and solving the optimality system consisting of the state system and the adjunct system that can be hard to solve. As an application, we considered a SIRS epidemic model with two controls whose expressions are given via selections of appropriately designed feedback map. These controls proved effective in reducing the number of infected individuals, either when used separately or when used together. Finally, as natural direction for our future works, we will propose extensions of our method to a class of epidemiological models where the dynamics are non-affine with respect to the control variable and also to control problems in the case of spatiotemporal epidemic models.

Footnotes

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

References

  • 1.Didier, R.: Les nouveaux risques infectieux: SRAS, grippe aviaire, et après?. Lignes de repères (2005)
  • 2.Bernoulli, D.: Essai d’une nouvelle analyse de la mortalité causée par la petite vérole et des avantages de l’inoculation pour la prévenir. Histoire de l’Acad. Roy. Sci.(Paris) avec Mém. des Math. et Phys. and Mém, pp. 1–45 (1760)
  • 3.Hethcote HW. The mathematics of infectious diseases. SIAM Rev. 2000;42(4):599–653. doi: 10.1137/S0036144500371907. [DOI] [Google Scholar]
  • 4.Keeling MJ, Rohani P. Modeling infectious diseases in humans and animals. Princeton University Press. 2008;42(4):599–653. [Google Scholar]
  • 5.Brauer F. Mathematical epidemiology: past, present, and future. Infect. Dis. Model. 2017;2(2):113–127. doi: 10.1016/j.idm.2017.02.001. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 6.Foppa IM. A Historical Introduction to Mathematical Modeling of Infectious Diseases: Seminal Papers in Epidemiology. Cambridge: Academic Press; 2016. [Google Scholar]
  • 7.Di Giamberardino P, Iacoviello D. Optimal control of SIR epidemic model with state dependent switching cost index. Biomed. Signal Process. Control. 2017;31:377–380. doi: 10.1016/j.bspc.2016.09.011. [DOI] [Google Scholar]
  • 8.Buonomo B, Lacitignola D, Vargas-De-León C. Qualitative analysis and optimal control of an epidemic model with vaccination and treatment. Math. Comput. Simul. 2014;100:88–102. doi: 10.1016/j.matcom.2013.11.005. [DOI] [Google Scholar]
  • 9.Zakary O, Rachik M, Elmouki I. On the analysis of a multi-regions discrete SIR epidemic model: an optimal control approach. Int. J. Dyn. Control. 2017;5(3):917–930. doi: 10.1007/s40435-016-0233-2. [DOI] [PMC free article] [PubMed] [Google Scholar]
  • 10.Mhlanga A. Dynamical analysis and control strategies in modelling Ebola virus disease. Adv. Differ. Equ. 2019;2019(1):458. doi: 10.1186/s13662-019-2392-x. [DOI] [Google Scholar]
  • 11.Moualeu DP, Weiser M, Ehrig R, Deuflhard P. Optimal control for a tuberculosis model with undetected cases in Cameroon. Commun Nonlinear Sci. 2015;20(3):986–1003. doi: 10.1016/j.cnsns.2014.06.037. [DOI] [Google Scholar]
  • 12.Xu D, Xu X, Xie Y, Yang C. Optimal control of an SIVRS epidemic spreading model with virus variation based on complex networks. Commun Nonlinear Sci. 2017;48:200–210. doi: 10.1016/j.cnsns.2016.12.025. [DOI] [Google Scholar]
  • 13.Bolzoni L, Bonacini E, Della Marca R, Groppi M. Optimal control of epidemic size and duration with limited resources. Math. Biosci. 2019;315:108232. doi: 10.1016/j.mbs.2019.108232. [DOI] [PubMed] [Google Scholar]
  • 14.Baba IA, Abdulkadir RA, Esmaili P. Analysis of tuberculosis model with saturated incidence rate and optimal control. Physica A. 2020;540:123237. doi: 10.1016/j.physa.2019.123237. [DOI] [Google Scholar]
  • 15.Abouelkheir I, El Kihal F, Rachik M, Elmouki I. Optimal impulse vaccination approach for an SIR control model with short-term immunity. Mathematics. 2019;7(5):420. doi: 10.3390/math7050420. [DOI] [Google Scholar]
  • 16.Berge T, Ouemba Tasse AJ, Tenkam HM, Lubuma J. Mathematical modeling of contact tracing as a control strategy of Ebola virus disease. Int. J. Biomath. 2018;11(07):1850093. doi: 10.1142/S1793524518500936. [DOI] [Google Scholar]
  • 17.Tan, J., Zou, X.: Optimal control strategy for abnormal innate immune response. Comput. Math. Methods Med. 2015 (2015) [DOI] [PMC free article] [PubMed]
  • 18.Gubar, E., Zhu, Q., Taynitskiy, V.: Optimal control of multi-strain epidemic processes in complex networks. In: International Conference on Game Theory for Networks, pp. 108–117 (2017)
  • 19.Omondi EO, Orwa TO, Nyabadza F. Application of optimal control to the onchocerciasis transmission model with treatment. Math. Biosci. 2018;297:43–57. doi: 10.1016/j.mbs.2017.11.009. [DOI] [PubMed] [Google Scholar]
  • 20.Pongsumpun P, Tang IM, Wongvanich N. Optimal control of the dengue dynamical transmission with vertical transmission. Adv. Differ. Equ. 2019;2019(1):176. doi: 10.1186/s13662-019-2120-6. [DOI] [Google Scholar]
  • 21.Pontryagin LS. Mathematical Theory of Optimal Processes. London: Routledge; 2018. [Google Scholar]
  • 22.Aubin JP, Bayen AM, Saint-Pierre P. Viability Theory: New Directions. Berlin: Springer; 2011. [Google Scholar]
  • 23.Kassara K. A unified set-valued approach to control immunotherapy. SIAM J. Control Optim. 2009;48(2):909–924. doi: 10.1137/07070591X. [DOI] [Google Scholar]
  • 24.Kassara K, Moustafid A. Angiogenesis inhibition and tumor-immune interactions with chemotherapy by a control set-valued method. Math. Biosci. 2011;231(2):135–143. doi: 10.1016/j.mbs.2011.02.010. [DOI] [PubMed] [Google Scholar]
  • 25.Gaff H, Schaefer E. Optimal control applied to vaccination and treatment strategies for various epidemiological models. Math. Biosci. Eng. 2009;6(3):469. doi: 10.3934/mbe.2009.6.469. [DOI] [PubMed] [Google Scholar]
  • 26.Kumar A, Srivastava PK. Vaccination and treatment as control interventions in an infectious disease model with their cost optimization. Commun. Nonlinear Sci. Numer. Simul. 2017;44:334–343. doi: 10.1016/j.cnsns.2016.08.005. [DOI] [Google Scholar]
  • 27.Boujallal L, Kassara K. State-input constrained asymptotic null-controllability by a set-valued approach. IET Control Theory Appl. 2015;9(15):2211–2221. doi: 10.1049/iet-cta.2014.1333. [DOI] [Google Scholar]
  • 28.Aubin JP. Viability Theory. Berlin: Springer; 2009. [Google Scholar]
  • 29.Aubin JP, Frankowska H. Set-Valued Analysis. Berlin: Springer; 2009. [Google Scholar]
  • 30.Deimling K. Multivalued Differential Equations. Berlin: Walter de Gruyter; 2011. [Google Scholar]
  • 31.Gutev V, Nedev S. Continuous selections and reflexive Banach spaces. Proc. Am. Math. Soc. 2001;129(6):1853–1860. doi: 10.1090/S0002-9939-00-05740-3. [DOI] [Google Scholar]
  • 32.Rachah, A., Torres, D. F.: Modeling, dynamics and optimal control of Ebola virus spread, vol. 1, no. 2, pp. 277–289 (2016). arXiv preprint arXiv:1603.05794
  • 33.Osemwinyen AC, Diakhaby A. Mathematical modelling of the transmission dynamics of ebola virus. Appl. Comput. Math. 2015;4(4):313–320. doi: 10.11648/j.acm.20150404.19. [DOI] [Google Scholar]

Articles from Journal of Applied Mathematics & Computing are provided here courtesy of Nature Publishing Group

RESOURCES