Epidemic management with admissible and robust invariant sets

We present a detailed set-based analysis of the well-known SIR and SEIR epidemic models subjected to hard caps on the proportion of infective individuals, and bounds on the allowable intervention strategies, such as social distancing, quarantining and vaccination. We describe the admissible and maximal robust positively invariant (MRPI) sets of these two models via the theory of barriers. We show how the sets may be used in the management of epidemics, for both perfect and imperfect/uncertain models, detailing how intervention strategies may be specified such that the hard infection cap is never breached, regardless of the basic reproduction number. The results are clarified with detailed examples.


Introduction
There is a large literature on the application of optimal control to epidemiology. Some of the earliest papers on the topic are [1], which investigates the optimal control of a disease of SIS type (susceptible-infective-susceptible), and [2], which considers an SIR (susceptible-infectiveremoved) model with vaccination to obtain optimal vaccination policies via dynamic programming. Other papers consider the optimal control of SIR models, [3][4][5][6]; of HIV [7,8]; and of malaria [9]. This list is by no means exhaustive, and we point the reader to the survey paper, [10]. Moreover, scores of papers have recently appeared involving the modelling and control of COVID-19, often via model-predictive control (MPC), see for example [11][12][13][14][15] and the thorough discussion of [16].
In this paper we build on the recent work concerning the application of set-based methods to epidemiology. We present a set-based analysis of two well-known compartmental epidemic models: the SIR and SEIR models, see for example [17,18], subjected to constraints on their inputs (social distancing, quarantine rate and/or vaccination rate) and state (a hard cap on the proportion of infectives), utilising the theory of barriers, [19,20]. These models often serve as good initial candidates to model new diseases, later being elaborated on to be more accurate.
It has only recently been shown that set-based methods may be used to maintain hard infection caps in epidemic models. To our knowledge, the paper [21] introduced the idea, describing the so-called viability kernel the set of all states for which there exists an input such that the a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 state and input constraints are satisfied for all future time, [22], of a two-dimensional model of a vector-borne disease. The work was then extended in [23] where the effects of modelling uncertainties of this system were investigated via the robust viability kernel. The recent paper [24] finds the viability kernel for an SIR model of a vector bourne disease by approximating the value function of an appropriate optimal control problem by solving the related Hamilton-Jacobi-Bellmann (HJB) equation.
In our previous research, [25], we showed that the theory of barriers may be utilised to describe the viability kernel (also known as the admissible set) as well as the maximal robust positively invariant set (MRPI) the set of states in which the state and input constraints are satisfied for all time regardless of the input of the malaria model considered in [21]. The theory of barriers describes the non-trivial part of these sets' boundaries, made of integral curves of the system that satisfy a minimum-/maximum-like principle, allowing them to easily be constructed. This contrasts with other approaches that estimate the viability kernel with algorithms that iteratively compute reachable sets, [26,Ch. 5], with interval analysis, [27], through the solution of HJB equations, [28], or via polynomial optimisation, [29], all of which suffer from the curse of dimensionality. Moreover, optimal control problems with hard state constraints (such as a hard cap on the number of infectives) are notoriously difficult to solve due to the presence of "active arcs", "jump conditions" and adjoint dynamics that are often difficult to specify. The construction of the sets via the theory of barriers does not suffer from these difficulties because the state constraints do not directly appear in their construction.
Other research that contributes to epidemiology with set-based ideas includes the paper [30], which describes the set of sustainable thresholds (in some sense the dual of the viability kernel) for the class of so-called cooperative epidemic models. An attractive aspect of this approach is that competing costs (for example, number of infective individuals versus the disease containment costs) may be analysed via Pareto efficiency. The work in [31,32] is closely related to set-based methods, but differs in that they do not describe or compute sets. Rather, these works combine Lyapunov-like functions with ideas from Viability theory (the so-called "regulation map") to specify continuous selections of set-valued maps that describe the elimination of a disease.
The set-based approach to the management of epidemics, as presented in [21,23,25], argues that intervention strategies (for example, the level of fumigation in the control of mosquito-borne diseases) should be chosen based on the location of the state with respect to the computed sets: If the state is located in the admissible set, then it is possible to maintain the infection cap for all time with a suitably chosen input (intervention strategy). If the state is located in the MRPI, then the cap will never be breached and any input may be used, such as, e.g., saving resources or minimising economic damage.
In the current paper we analyse these sets for the constrained SIR and SEIR models under two cases: a perfect model case and an imperfect model case. For the perfect model case the modelling parameters are assumed to be perfectly known and the intervention strategy may be chosen depending on the location of the state, as was described earlier. For the imperfect model case we assume that there is a pre-designed intervention strategy, possibly a feedback, that has been designed and implemented on the system. For this case we only describe the MRPI, which corresponds to states from which the infection cap can always be maintained by the intervention strategy regardless of the modelling uncertainty. This latter aspect could be very valuable to epidemiologists because the parameters appearing in epidemic models sometimes cannot be accurately estimated (see e.g. [16]).
The contributions of the paper are as follows: • To our knowledge, this is the first paper that deals with the problem of maintaining a hard infection cap for the SIR and SEIR epidemic models (under the assumption of both perfect and imperfect modelling), via a set-theoretic approach: the theory of barriers, [19,20].
• With the theory of barriers we are not only able to decrease the dimension in the computation of the sets, but also able to obtain easily checkable inequalities of the systems' parameters that allow one to classify the two sets, for the considered models, into trivial or nontrivial ones.
• We demonstrate that the MRPI, whose use in epidemics was introduced in our previous research [25], is also useful as a tool to maintain an infection cap for uncertain epidemic models.
• A note-worthy observation made in the paper regards the basic reproduction number, R 0 . In epidemiology intervention strategies are often designed so that R 0 is less than one, which guarantees that the disease will eventually die out. We demonstrate that using the computed sets the infection cap can be maintained regardless of the value of R 0 , even if it is greater than one.
The paper is organised as follows. We present the SIR and SEIR models we will study in Section 2, along with precise formulations of the three problems we will address in the paper. We summarise the relevant results from the theory of barriers in Section 3, which we apply to describe the sets for the two epidemic models in Sections 4, 5 and 6. In Section 7 we present conditions involving the system parameters that allow one to easily determine whether the sets are trivial or not. Detailed examples are presented in Section 8 and the paper is concluded with a discussion in Section 9.

Two epidemic models
In this section we describe the two well-known epidemic models on which the paper will focus: the SIR and SEIR models.

SIR model
The (normalised) SIR model, see for example [17,Ch. 2] and [18], is described by the following ordinary differential equations: We assume that the birth and death rates are zero. In other words, we assume that the disease is relatively short-lived such that the effects of demographic change can be ignored. Therefore, assuming an initial state satisfying S(t 0 ) + I(t 0 ) + R(t 0 ) = 1, all integral curves of the system remain in the set: for all t � t 0 . Let N 2 R >0 denote the constant population size, and let s ≜ NS, i ≜ NI and r ≜ NR denote the number of susceptibles, infectives and recovered individuals, respectively. Let b 2 R >0 , the contact rate, denote the average number of contacts a person makes per unit time that is sufficient for disease transmission (this contact can be with infectives and/or other susceptibles). Then βI is the average number of contacts with infectives a susceptible makes per unit time, and βIs is the number of new infectives per unit time. Thus, βIS is the proportion of new infectives per unit time. It is assumed that recoveries occur according to a Poisson process with arrival rate g 2 R >0 . See [17,Ch. 2], [18] and [13] for a more in-depth discussion of this model's derivation.
In conclusion, β may be interpreted as the average number of people a member of the population comes into contact with over one day, and 1 g is the average number of days an infective individual remains infective. [0, 1] denoting exposed individuals. This is a general model for diseases where individuals only become infective after a non-negligible period of time from exposure. It is assumed that an individual's evolution from being exposed to being infective follows a Poisson process with arrival rate Z 2 R >0 . Thus, 1 Z is the average number of days it takes an exposed individual to become infective.

SEIR model
As before, we assume that the total population remains constant and that the SEIR model is normalised. Thus, all solutions remain in the (redefined) set for all t 2 [t 0 , 1[.

On set-based methods applied to epidemiology
Now we present an introductory overview of how we intend to apply set-based methods to the management of epidemics. The subsequent sections will present rigorous treatment of the ideas. Consider the SIR model, (1)-(3), and impose a hard infection cap, I(t) � I max , for all time. Moreover, assume that the contact rate, β, may be changed over time within some bounds, representing social distancing measures. That is, assume β(t)2[β min , β max ], 0 < β min � β max where β min is the minimal contact rate (and thus represents the maximal social distancing) and β max is the maximal or nominal contact rate (and may thus represent minimal social distancing). Consider the SEIR model, (5)- (8), and also impose a hard infection cap, I(t) � I max , for all time. Assume that in addition to the contact rate, β, the parameter γ may also be changed over time within some bounds, i.e., γ(t) 2 [γ min , γ max ], 0 < γ min � γ max . This may represent quarantine or vaccination measures.
The first problem we want to address is the following: Problem 1. For the SIR model, supposing the parameter γ is perfectly known (SIR prefect model), how should the contact rate β(t), be specified over time such that the infection cap, I max , is never breached? For the SEIR model, supposing the parameter η is perfectly known (SEIR prefect model), how should the contact rate β(t) and the removal rate γ(t) be specified over time such that the infection cap is never breached?
Our approach to solving this problem will be to find the admissible set, A, which, for the SIR model, will correspond to all states from which there exists an input β satisfying β(t)2 [β min , β max ] for all t 2 [t 0 , 1[ and such that the infection cap, I max , is never breached. For the SEIR model, A will correspond to all states for which inputs β(t) 2 [β min , β max ] and γ(t) 2 [γ min , γ max ] exist such that the infection cap is never breached. Then the contact rate β, as well as the removal rate γ for the SEIR model, will be deduced depending on the location of the state with respect to the computed sets. Moreover, the complement of the admissible set, A C , is the set of states for which the infection cap cannot be maintained, regardless of the interventions. Thus, from this set the problem has no solution and if the state is located here, either the infection cap needs to be relaxed, or the maximal allowable social distancing measures (β min ), and/or maximal quarantine rate for the SEIR model (γ max ), need to be strengthened.
The set A is obtained by constructing its boundary, which will be seen to be made of two parts, the usable part, a subset included in the boundary of the state constraints, and the barrier, entirely contained in the interior of the state constraint set (see the general presentation of Section 4, and further results and comparisons for the SIR and SEIR perfect models in Sections 7 and 8).
The second problem concerns also the SIR and SEIR perfect models: These sets are called the Maximal Robust Positively Invariant (MRPI) sets associated to the SIR and SEIR perfect model respectively (see their general definitions 3.2 and 3.3). Their construction, again via the construction of their boundary, consisting, as before, of a usable part and a barrier, is presented in Sections 5, 7 and 8).
The third problem concerns the SIR and SEIR imperfect models: Problem 3. Suppose now that, for the SIR model, the value of the parameter γ is not known.
All that we know is that γ 2 [γ min , γ max ]. We also assume that an intervention strategybðtÞ : ½t 0 ; 1½! ½b min ; b max � has been pre-designed (SIR imperfect model). Find the states from which this intervention strategy will maintain the infection cap regardless of the error in the parameter γ, for all t 2 [t 0 , 1[. For the SEIR model, we pose the same problem with unknown parameter η 2 [η min , η max ] in place of γ and with predesigned intervention strategiesbðtÞ : ½t 0 ; 1½! ½b min ; b max � andĝðtÞ : ½t 0 ; 1½! ½g min ; g max � (SEIR imperfect model).
Its solution will consist in constructing the MRPI set associated to the SIR (resp. SEIR) imperfect model (see Sections 6,7 and 8). Table 1 summarises the status of the parameters and controls for each problem. The main mathematical tool to solve these three problems is sketched in the next section.

The theory of barriers
We now summarise the relevant theory from [19,20], which we will apply to describe the sets of the two epidemic models. Consider a state and input constrained nonlinear system: where xðtÞ 2 R n is the state; x 0 is the initial condition; uðtÞ 2 R m is the input; U is the set of Lebesgue measurable functions that map the interval ½t 0 ; 1½� R to a compact and convex set U � R m (note the difference between U, a function space, and U, a subset of R m ); and g i , i = 1, . . ., p, are the state constraints. We impose the following assumptions: (A2) Every x 0 , with kx 0 k<1, and every u 2 U admits a unique absolutely continuous integral curve of (10) that remains bounded over any finite interval of time.
(A4) The function g i , i = 1, . . ., p, is C 2 and the set of points given by g i (x) = 0 defines an n − 1 dimensional manifold.
It can be verified that the SIR and SEIR models analysed in subsequent sections satisfy Assumptions (A1)-(A4) for each considered interpretation of the input. In particular, the following growth condition is always satisfied: there exists a C < 1 such that sup u2U jx > f ðx; uÞj � Cð1þ k x k 2 Þ for all x 2 R n , which implies (A2). An important consequence of these assumptions is that the admissible and MRPI sets are closed. We use the following notation: let x ðu;x 0 ;t 0 Þ denote the unique solution of (10) from the initial state x 0 2 R n at initial time t 0 , with input function u 2 U, and let x ðu;x 0 ;t 0 Þ ðtÞ denote the solution at time t 2 [t 0 , 1[. If clear from context, we will use the notation x u and x u (t) instead. Let G ≜ fx : g i ðxÞ � 0; i ¼ 1; 2; :::; pg, G À ≜ fx : g i ðxÞ < 0; i ¼ 1; 2; :::; pg, and G 0 ≜ fx : 9i 2 f1; 2; :::; pg s:t: g i ðxÞ ¼ 0g: The set of indices of active constraints at a point x is denoted IðxÞ ¼ fi : g i ðxÞ ¼ 0g. We let L f gðx; uÞ ≜ DgðxÞf ðx; uÞ denote the Lie derivative of a differentiable function g : R n ! R with respect to f(�, u) at the point x.
We now present the definitions of the two sets that we will describe for the compartmental epidemic models.
Definition 3.1. The admissible set also known as the viability kernel in Viability theory, [22], see [19], of the system (10) and (11), denoted by A, is the set of initial states for which there exists a u 2 U such that the corresponding solution satisfies the constraints (11) for all future A ≜ fx 0 2 G : 9u 2 U; x ðu;x 0 ;t 0 Þ ðtÞ 2 G 8t 2 ½t 0 ; 1½ g:  (10) and (11) contained in G, see [20], is the union of all RPIs that are subsets of G. Equivalently this equivalence is shown in [20], We obviously have M � A. Moreover, under (A1)-(A4) the sets A and M are closed, and we denote their boundaries by @A and @M, respectively.
The main result from the papers [19,20] is a characterisation of the sets' boundaries that are made up of two complementary parts: the usable part, which is contained in the boundary of the constrained state space, ½@A� 0 ≜ @A \ G 0 and ½@M� 0 ≜ @M \ G 0 ; and the barrier, which is contained in the interior of the constrained state space, ½@A� À ≜ @A \ G À and ½@M� À ≜ @M \ G À . Points on the admissible set's usable part satisfy: whereas those on the MRPI's usable part satisfy: The barrier is made up of integral curves of the system that satisfy a minimum-/maximum-like principle that may intersect the boundary of the constraint set, G 0 , in a tangential manner. We summarise these results in the following theorem, for which proofs can be found in [19,20].
Theorem 3.1. Under the assumptions (A1)-(A4), every integral curve x � u contained in ½@A� À (resp. ½@M� À ) and crossing G 0 tangentially in finite time satisfy, together with the corresponding control � u 2 U, the following necessary conditions. There exists a nonzero absolutely continuous maximal solution l � u to the adjoint equation: such that the Hamiltonian λ T f satisfies: for almost all t. Moreover, we have: To now compute the curves running along the barrier, we first obtain the points of ultimate tangentiality via the expressions (17) (resp. (18)) and then integrate the system backwards in time with the control function that minimizes (in the case of the set A) or maximizes (in the case of the MRPI set M) the Hamiltonian for almost every t, according to the expression (15) (resp. (16)) with λ satisfying (14).
In the sequel the set G will be expressed by a single state constraint that caps the proportion of infective individuals, and we will interpret the function u in two ways. The first interpretation will be that u is a time-varying "controllable" input (u = β for the SIR model and u = (β, γ) for the SEIR model), which will allow us to solve Problem 1 (resp. Problem 2) (see Subsection 2.3 and Table 1), i.e., find the contact rate β (in the SEIR case also the removal rate γ) such that, the parameter γ being perfectly known, the infection cap, I max , is never breached (resp. never breached for all β 2 [β min , β max ]).
The second interpretation will be that u is an "uncontrollable" disturbance term (u = γ for the SIR model and u = η for the SEIR model) encapsulating uncertainty in modelling parameters. This will allow us to address Problem 3, i.e., given a pre-designed feedback, find the states from which the infection cap will be always maintained regardless of the error in the modelling parameters.

Carrying out the analysis of the admissible set for the perfect SIR and SEIR models (problem 1)
In this section we carry out the detailed analysis of the admissible sets of the constrained SIR and SEIR models under the assumption of perfect modelling (Problem 1). For the sake of readability, comparisons of admissible and robust sets for the perfect and imperfect models as well as numerical applications are presented in Sections 7 and 8.

Admissible set for perfect SIR model
We assume that the SIR model is perfect (that is, γ is known perfectly) and that β is a controllable time-varying input, which may, for example, model social-distancing measures to halt the disease's spread. We impose hard constraints on the input and state: 1] is the hard infection cap, and b min ; b max 2 R >0 , with 0 < β min � β max , are the minimal and maximal contact rates, respectively. Recall that the bound β max may be interpreted as the nominal contact rate in absence of the disease, and β min as the maximal social distancing.
Because R has no influence on S or I we can study the two-dimensional subsystem: _ IðtÞ ¼ bðtÞSðtÞIðtÞ À gIðtÞ; ð20Þ with t 2 [t 0 , 1[. To ease our notation we will sometimes label the state x ≜ ðS; IÞ > 2 R 2 , the control u ≜ b 2 R, the control constraint set U ≜ ½b min ; b max � � R, the single constraint function gðxÞ ≜ I À I max , and the set G P ≜ G \ P, where P is given by (4). We also denote by Applying condition (12), the usable part ½@A� 0 is given by the condition Therefore the usable part is the set Focusing now on the barrier ½@A� À , let us label tangent points to G 0 by where � t denotes the time at which G 0 is reached, as explained in Theorem 3.1. Invoking ultimate tangentiality, (17), we get: Dg(z 1 , z 2 ) = (0, 1) and ð22Þ and along with the constraints z 2 = I max and z 1 + z 2 � 1, we get that the single tangent point must satisfy (β min z 1 − γ)I max = 0, or ðz 1 ; z 2 Þ ¼ ð g b min ; I max Þ, with g b min þ I max � 1, one of the end point of the usable part.
Proposition 4.1. The input � b, given by (25), associated with the barrier of the admissible set of the perfect SIR model, (19)- (21) only switches at isolated points in time.
The next Proposition states that barrier curves associated with the admissible set of the perfect SIR model always evolve backwards into G − . Proposition 4.2. Consider the constrained perfect SIR model, (19)- (21). The barrier ½@A� À is constituted by the curve x � u generated by (25) with (23), that ends tangentially to G 0 \ P at the point ð g b min ; I max Þ 2 G 0 \ P with g b min þ I max � 1, and is contained in G − at least for all t 2� � t À �; � t½, � > 0.
Proof. Let � g denote the mapping t ! � g ðtÞ ≜ gðx � u ðtÞÞ ¼ IðtÞ À I max over some open interval � � t À �; � t½, with � > 0. We have shown that the barrier curve x � u ¼ ðSðtÞ; IðtÞÞ has to be such that � g vanishes at � t as well as its first derivative, which is equal to d� g dt ¼ L f g. Moreover, β must remain constant, equal to β min , in the interval � � t À �; � t½. Let us compute � g 's left second derivative: Using the fact that bð � tÞSð � tÞIð � tÞ À gIð � tÞ ¼ 0, Sð � tÞ ¼ g b min , bð � tÞ ¼ b min and Ið � tÞ ¼ I max , this simplifies to: Therefore, using Taylor's expansion formula, in an interval � � t À �; � t½ for a sufficiently small �, which proves that the barrier curve x � u remains in G − in this interval.
Remark 4.2. Note that the point of the barrier satisfying I = 0 is an equilibrium point: Therefore the unique solution from (S(t 0 ), 0) is reduced to this point.
In fact, the whole axis I = 0 may be seen to be included in ½@A� À since every point of this axis is an equilibrium and thus satisfies the constraint forever and belongs to @P, the boundary of P (see also Section 7.1). It is remarkable that the points in this axis, though in the barrier, do not satisfy the conditions of Theorem 3.1 since they do not cross G 0 (recall that Theorem 3.1 is only necessary).
Finally, the integral curve of the system S(t) = 0, IðtÞ ¼ e À gð� t À tÞ I max , t 2� À 1; � t�, i.e. the semiopen line segment {(0, I)jI 2 ]0, I max ]} also lies in A \ @P and hence in ½@A� À , neither satisfies the conditions of Theorem 3.1 since it does not intersect G 0 tangentially. It results that ½@A� À is the union of the set made by the integral curves satisfying Theorem 3.1 and of a set of special curves that do not intersect G 0 tangentially.
Suppose now that λ 3 (t) = 0 over some nonempty interval I. Then we would have _ l 3 ðtÞ ¼ À bðtÞSðtÞðl 2 ðtÞ À l 1 ðtÞÞ ¼ 0 over this interval, which would imply that S(t)(λ 2 (t) − λ 1 (t)) = 0, provided that S(t)6 ¼0, over this interval. But if S(t) = 0 over some interval of time, according to the uniqueness of the solution of (26)- (29), this would imply that S(t) = 0 for all times and in particular Sð � tÞ ¼ z 1 ¼ 0, which is excluded by assumption. Therefore λ 3 (t) = 0 over some nonempty interval would imply λ 2 (t) − λ 1 (t) = 0 over this interval, which we have established is impossible in the first part of the proof.
The next proposition identifies which parts of T A are associated with barrier curves that evolve backwards into G − . To now numerically construct the admissible set one only needs to focus on these tangent points.

Carrying out the analysis of the MRPI for the perfect SIR and SEIR models (problem 2)
In this section we carry out the detailed analysis of the Maximal Robust Positively Invariant sets of the constrained SIR and SEIR models under the assumption of perfect modelling (Problem 2).
As for the previous section, comparisons and numerical applications are presented in Sections 7 and 8.

MRPI for perfect SIR model
We consider the same setting as that of Subsection 4.1, the constrained two-dimensional system (19)- (21), and describe the system's MRPI set.
We first remark that the usable part ½@M� 0 , according to condition (13), is given by max b2½b min ;b max � L f g ¼ bSI À gI � 0; ðS; IÞ 2 G P s:t: I À I max ¼ 0; Now, carrying out the analysis with the ultimate tangentiality condition, (18), and Hamiltonian maximisation condition (16), we note that the only difference with respect to (22) is that we maximize with respect to β in place of minimizing. Thus, we find that the sole tangent point is located at The adjoint equation being the same as (23), the input associated with the MRPI barrier curve is given by: anything if l 2 ðtÞ À l 1 ðtÞ ¼ 0; which also only switches at isolated points in time, the proof being exactly the same as that of Proposition 4.1.
The following proposition is the analogue of Proposition 4.2. Its proof follows exactly the same lines and is left to the reader.  (21). The barrier ½@M� À is constituted by the curve x � u generated by (34) with (23), that ends tangentially to G 0 \ P at the point ð g b max ; I max Þ 2 G 0 \ P with g b max þ I max � 1, and is contained in G − at least for all t 2� � t À �; � t½, � > 0.
The reader may also verify that Remark 4.2 also applies to the barrier of the MRPI set.

MRPI for perfect SEIR model
We consider the same setting as that of Subsection 4.2: the constrained three-dimensional system (26)-(29) and describe the system's MRPI set. The analysis follows exactly the same lines as in Subsection 4.2, the operation of maximisation replacing the minimisation one. We therefore only sketch the results. The usable part ½@M� 0 is given by max g2½g min ;g max � L f g ¼ ZE À gI � 0; ðS; E; IÞ 2 G P s:t: I À I max ¼ 0: Therefore, Concerning the barrier, using conditions (16) and (18), we identify the line segment of potential tangent points: where the subscript M P in T M P refers to the "perfect" model case.
The adjoint system being described by (31), the input associated with the MPRI's barrier curves are thus given by: anything if l 3 ðtÞ ¼ 0: These inputs may only switch at isolated points in time provided that z 1 6 ¼ 0 in (35), the proof following exactly the same arguments as the proof of Proposition 4.3. The next proposition is the analogue of Proposition 4.4, its proof also following the same lines of reasoning.
Proposition 5.2. Consider the constrained perfect SEIR model, (26)- (29). The barrier ½@M� À is constituted by the integral curves x � u satisfying (31) and (36), that end tangentially to G 0 at a point of T M P , defined by (35), for which and is contained in G − at least for all t 2� � t À �; � t½, � > 0. As before, the reader may verify that Remark 4.3 applies to the MPRI for the perfect SEIR model as well.

MRPI for the imperfect SIR and SEIR models (problem 3)
In this section we carry out the detailed analysis of the MRPIs of the constrained SIR and SEIR models under the assumption of imperfect modelling, corresponding to Problem 3).

MRPI for imperfect SIR model
We again consider the normalised SIR model, but now assume that the model parameter γ is not perfectly known and that the intervention strategy, a simple affine feedback strategy, is given by: where the level of social distancing depends on the proportion of infective individuals, has been pre-designed and implemented. Thus, the system is now: _ IðtÞ ¼bðtÞSðtÞIðtÞ À gðtÞIðtÞ; IðtÞ À I max � 0; gðtÞ 2 ½g min ; g max �; ð41Þ with t 2 [t 0 , 1), and 0 < γ min � γ max . Again, to ease our notation, we will sometimes label the state x ≜ ðS; IÞ > 2 R 2 , the uncertain model parameter u ≜ g 2 R, the parameter set U ≜ ½g min ; g max � � R, and gðxÞ ≜ I À I max .
The adjoint Eq (14) reads: where aðIÞ ≜ 2 b min À b max Invoking conditions (16) we find the input associated with the MRPI barrier curve, given by: anything if l 2 ðtÞ ¼ 0: We again have a result concerning singular barrier curves. Proposition 6.1. Consider the constrained imperfect SIR model, (39)-(41). Along any barrier integral curve satisfying (42) and (43) and such that z 1 6 ¼ 0 and 2β min < β max , the input � g, given by (43), only switches at isolated points in time.
From condition (18) the sole tangent point is located at z ¼ ð g min b min ; I max Þ, with g min b min þ I max � 1, which we derive from the fact thatbðI max Þ ¼ b min .
Again, the barrier curve evolves backwards into G − : Proposition 6.2. Consider the constrained imperfect SIR model, (39)-(41). The barrier ½@M� À is constituted by the curve x � u generated by (43) with (42), that ends tangentially to G 0 \ P at the point ð g min b max ; I max Þ 2 G 0 \ P, with g min b min þ I max � 1, and is contained in G − at least for all t 2� � t À �; � t½, � > 0.

MRPI for imperfect SEIR model
Similar to the previous case, we consider the system: _ EðtÞ ¼bðtÞSðtÞIðtÞ À ZðtÞEðtÞ; ð45Þ IðtÞ À I max � 0; ZðtÞ 2 ½Z min ; Z max �; ð47Þ with t 2 [t 0 , 1[, 0 < η min � η max . The inputsb : ½t 0 ; 1½! ½b min ; b max � andĝ : ½t 0 ; 1½! ½g min ; g max � are assumed to be pre-designed by the simple affine feedback strategies: andĝ which may be interpreted as the fact that we increase social distancing and increase the rate of quarantining if the infective number increases. For this setting we identify: x ≜ ðS; E; IÞ > 2 R 3 , the uncertain modelling parameter u ≜ Z 2 R, the uncertain parameter set U ≜ ½Z min ; Z max � � R, gðxÞ ≜ I À I max , and G P ≜ G \ P.
As before, the reader may easily verify that the set of potential tangent points is: where the subscript M I in T M I denotes the "imperfect" case. The adjoint, (14), reads: where dðIÞ ≜ 2 g max À g min I max � � I þ g min , and α(I) is as in Subsection 6.1.
It can be verified that the input associated with barrier curves of the MRPI is given by: Z min if l 3 ðtÞ À l 2 ðtÞ < 0; anything if l 3 ðtÞ À l 2 ðtÞ ¼ 0: , the input � Z, given by (52), only switches at isolated points in time. Proof. Assume, as before, that λ 3 (t) − λ 2 (t) = 0 in a given interval of time I. We immediately deduce that _ Thus, if the above determinant does not vanish for all I 2 [0, I max ] and if S 6 ¼ 0 we arrive at λ 1 = λ 2 = λ 3 = 0, the desired contradiction.
The proof of the next proposition follows the same lines as before and is left to the reader. Proposition 6.4. Consider the constrained imperfect SEIR model, (44)-(47). The barrier ½@M� À is constituted by the curves x � u generated by (52) with (51), that end tangentially to G 0 \ P at a point of T M I , defined by (50), with and is contained in G − at least for all t 2� � t À �; � t½, � > 0.

Further results on the relative location of the admissible and MRPI sets
The admissible and/or MRPI sets may be trivially equal to the entire constrained state space, G P , in which case any intervention strategy will maintain the infection cap. We now present some results that summarise the relationship between certain inequalities of the system parameters and the relative location of the sets.  (21). We have the following:

SIR model
Consider the constrained imperfect SIR model (39)-(41), under the feedback (38). We have the following: Proof. The constraint set G P is the polyhedron: and its boundary is the union of the four faces The last three inequalities (55)-(57) are trivially satisfied and (54) results from (53), which proves that (53) implies M ¼ G P and thus also M ¼ A ¼ G P . The converse, namely that M ¼ G P implies (53), follows from the fact that (54) is valid only if (53) holds.
We proceed in the same way to prove that A ¼ G P if b max > g ð1À I max Þ and b min � g ð1À I max Þ , up to the fact that, in (54)-(57), the maximisation with respect to β is replaced by the minimisation with respect to β.
The third bullet point results from the fact that none of the two previous cases hold, which prevents from having an equality of M or A to G P .
The last two bullet points are obtained as the first and third bullets by changing f in f ðS; I; gÞ ¼ ÀbSÎ bSI À gI ! and maximising with respect to γ 2 [γ min , γ max ].
Proposition 7.1 provides easily checkable inequalities of the system parameters that may be used to determine whether the sets are trivial or not. For example, focusing on the perfect SIR model, if b max � g ð1À I max Þ then G P is invariant for any choice of input and one need not worry about ever breaching the infection cap, as long as the state initiates inside G P . However, this is not the case when b max > g ð1À I max Þ , and the barrier of the MRPI must be computed by integrating the coupled differential equations of the system and the adjoint ones. If b min > g ð1À I max Þ then the MRPI and the admissible set are proper subsets of G P . In this case there are parts of G P from which the infection cap will definitely be breached and, starting from A or M, we must compute the corresponding barrier to derive an adapted intervention strategy.

SEIR model
We now summarise the analogous results for the perfect and imperfect SEIR models. Proposition 7.2. Consider the constrained perfect SEIR model (26)- (29). We have the following: Consider the constrained imperfect SEIR model (44)-(47) under the feedback (48) and (49). We have the following: Proof. We follow the same arguments as in the proof of Proposition 7.1. For the SEIR model (26)- (29), the set G P is the polyhedron and @G P ¼ fðS; E; IÞ 2 G P j I ¼ I max or S ¼ 0 or E ¼ 0 or I ¼ 0 or S þ E þ I ¼ 1g: We thus prove that M ¼ G P is equivalent to the first bullet point inequality, i.e.
Therefore, f points towards the interior of G P along @G P for all values of β and γ, or equivalently M ¼ A ¼ G P , if and only if (58) holds.
We proceed in the same way to prove that A ¼ G P if and only if η(1 − I max ) − γ min I max > 0 and η(1 − I max ) − γ max I max � 0, up to the fact that the maximisation with respect to β and γ is replaced by the minimisation with respect to β and γ.
As in the proof of Proposition 7.1, the third bullet point results from the fact that none of the two previous cases hold, which prevents from having an equality of M or A to G P .
Finally, the last two bullet points are obtained as the first and third bullets by changing f in  (21), with γ = 0.5, β min = 0.6, β max = 0.8 and I max 2 {0.02, 0.15, 0.4}. Thus, this (hypothetical) disease has an average recovery rate of 2 days and, on average, an individual nominally comes into contact (sufficiently to catch the disease) with 0.8 other people per day. We have chosen these parameters to demonstrate the different types of possible sets.
For I max = 0.4 we see from Proposition 7.1 that M ¼ A ¼ G P . Thus, with such a large cap on the infection numbers the entire space G P is robustly invariant. For I max = 0.15 we see from Proposition 7.1 that M ⊊ A ¼ G P . To now obtain the barrier associated with the MRPI, we integrate the system backwards in time from the tangent point z ¼ ð g b max ; I max Þ utilising the input (34). We stop integrating once the trajectory intersects the boundary of G P because the set P is positively invariant (recall the discussions in Subsections 2.1 and 2.2.) For I max = 0.02 we see from Proposition 7.1 that M ⊊ A ⊊ G P . To obtain the barriers of the admissible set and MRPI we integrate the system backwards in time from the tangent points z ¼ ð g b max ; I max Þ and z ¼ ð g b min ; I max Þ utilising the input (34) and (25), respectively. See Fig 1 for details. It turns out that the input associated with the admissible set is saturated at β min , and the one associated with the MRPI is saturated at β max . Moreover, we can see that, going backwards, the barrier curves always evolve into G − , as is expected from Propositions 4.2 and 5.1.
To now use the sets in the management of an epidemic the intervention strategy may be chosen according to the location of the state. A possible strategy could be: • If ðSðtÞ; IðtÞÞ > 2 @A, let β(t) = β min , • If ðSðtÞ; IðtÞÞ > 2 ½@A� 0 , let bðtÞ ¼ g SðtÞ .
To clarify, if the state is located in the MRPI then it is guaranteed that the infection cap can always be maintained, and the nominal contact rate, β max , may be allowed. The same is true if the state is located in the interior of the admissible set, but this may only be possible for some period of time. If the state reaches the admissible set's barrier, ½@A� À , then maximal social distancing must be imposed, thus β(t) = β min , and the infection cap will be reached in finite time, but never breached. If the state reaches the admissible set's usable part, ½@A� 0 , then some freedom in the social distancing is allowed. In fact, we may choose β(t) such that _ I ¼ 0, which yields bðtÞ ¼ g SðtÞ . Lastly, the set A C should always be avoided, as from here the infection cap is guaranteed to be breached regardless of β. Fig 2 shows (29), with the constants β min = 0.8, β max = 1, g min ¼ 1 5 , g max ¼ 1 3 , Z ¼ 1 5 and I max 2 {0.3, 0.4}. From Proposition 7.2 we see that M ⊊ A ¼ G P for I max = 0.4. Thus, we sample final tangent points along the line segment T M P for which z 1 � minf g min b max ; 1 À g min I max =Z À I max g (see Proposition  For I max = 0.3 we see that M ⊊ A ⊊ G P . Thus, for this case we find both M and A, integrating backwards from T M P and T A with the appropriate inputs as given in Subsections 5.2 and 4.2. The inputs � b and � g are saturated at β max and γ min (resp. β min and γ max ) for the barrier of the MRPI (resp. the barrier of the admissible set). As was observed in the examples for the SIR model, if I max is large enough the entire set G P is robustly invariant, and if I max is small enough, the sets A and M become nontrivial. 8.2.2 Imperfect SEIR example. Now consider the imperfect SEIR model, (44)-(47) with the parameters β min = 0.8, β max = 1, g min ¼ 1 5 , g max ¼ 1 3 , Z min ¼ 1 7 , Z max ¼ 1 5 and I max = 0.1, and the feedback (48) and (49). From Proposition 7.2, we have M ⊊ G P . Thus, we integrate backwards from points in T M I for which z 1 � minf g max b min ; 1 À g max I max =Z max À I max g (see Proposition 6.4) with the input � Z as specified in Subsection 6.2. An unexpected result for this set is that the input associated with barrier curves, � Z, switches from η min to η max at some point along the barrier (these points are indicated on the figure). To explain this, we can interpret η as an input with the goal of maximising I over time. For this example it is optimal for it to remain at η min in order to build up the proportion of exposed individuals to some threshold when it switches to η max , which then results in a large increase in the proportion of infectives. This is the worst-case change in η that can be expected. The computed set is shown in Fig 5.

Discussion
We applied the theory of barriers to describe integral curves of the system that run along special parts of the boundaries of the admissible and MRPI sets of the constrained SIR and SEIR models. The analysis in Sections 4, 5 and 6 summarises the details required to construct the sets for each of the three problems stated in Subsection 2.3, with Section 7 presenting inequalities of the system parameters that allow one to determine whether the sets are trivially equal to the constrained state space G P or not. We also demonstrated, in the example in Subsection 8.1 how an intervention strategy may be chosen using the computed sets in order to maintain the infection cap.
There are many avenues of future research to pursue: the analysis could be applied to other epidemic models; we could investigate the efficacy of combining model predictive control with knowledge of the sets; or look into the effects of more realistic intervention strategies that must remain constant over long periods of time, as opposed to the continuous ones introduced in (38) or (49)-(48).