key: cord-0930476-2xflheth
authors: Yan, Xiefei; Zou, Yun
title: Optimal and sub-optimal quarantine and isolation control in SARS epidemics
date: 2007-05-18
journal: Math Comput Model
DOI: 10.1016/j.mcm.2007.04.003
sha: f91e46f590c6f6c132760ffd3414c6b92119e093
doc_id: 930476
cord_uid: 2xflheth

This paper discusses the application of optimal and sub-optimal controls to a SEQIJR SARS model via the Pontryagin’s Maximum Principle. To this end, two control variables representing the quarantine and isolation strategies are considered in the model. The numerical optimal control laws are implemented in an iterative method, and the sub-optimal solution is computed using a genetic algorithm. The simulation results demonstrate that the maximal applications of quarantining and isolation strategies in the early stage of the epidemic are of very critical impacts in both cases of optimal and sub-optimal control. Otherwise, the control effect will be much worse. This gives a theoretical interpretation to the practical experiences that the early quarantine and isolation strategies are critically important to control the outbreaks of epidemics. Furthermore, our results also show that the proposed sub-optimal control can lead to performances close to the optimal control, but with much simpler strategies for long periods of time in practical use.

Severe acute respiratory syndrome (SARS), a new, highly contagious, viral disease, emerged in China late in 2002 and quickly spread to 32 countries and regions, causing in excess of 774 deaths and 8098 infections worldwide [1] . Since there are no valid medicines or vaccines for SARS, measures to control the spread of SARS had to take two major forms: isolation of symptomatic individuals and quarantining and close observation of asymptomatic individuals [2] . It had been shown that quarantining and isolating diseased individuals is a critically important strategy that can control SARS outbreaks, because of a reduction in the contact rate between susceptible and diseased individuals [2] [3] [4] . Moreover, Gumel [3] examines the impact of isolation and quarantine on the control of SARS during the outbreaks in Toronto, Hong Kong, Singapore and Beijing in a mathematical way. It uses a deterministic model that closely mimics the data for cumulative infected cases and SARS-related deaths.

In addition, mathematical models have recently been used to examine the transmission dynamics and to model the control of SARS in the literature [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] . But most of them are mainly in the purely simulative levels [5] [6] [7] [8] [9] [10] .

Wang and Ruan [9] propose a mathematical model consists of six subpopulations, namely susceptible, exposed, quarantined, suspect, probable, and removed, to simulate the SARS outbreak in Beijing. Chowell and Castillo-Chavez [10] use the uncertainty and sensitivity analysis of the basic reproductive number R 0 to assess the role that the model parameters play in outbreak control, and the results show that the transmission rate and isolation effectiveness have the largest fractional effect on R 0 . Shi [11] and Lloyd-Smith [12] present a stochastic dynamic model of SARS spreading, and a stochastic model of a SARS outbreak in a community and its hospital, respectively. Furthermore, small-world networks have been introduced to study the propagation of SARS by Naoki Masuda [13] and Small [14] . It has been shown that if all infectious individuals are isolated as rapidly as they are identified, the severity of the outbreak would be minimal [14] . Moreover, Gumel [3] and Chowell [4] attempt to obtain a threshold for the basic reproductive number R 0 for assessing the strategies of quarantine and isolation, and discuss the control of SARS by looking at the role of disease transmission parameters in the reduction of R 0 and the prevalence of the disease, while the discrete case has been studied by Zhou and Ma [15] . However, those models do not consider control strategies, since their discussions are based on prevalence of disease at equilibria. In Swan's work [16] , optimal control theory is applied to obtain maximal benefits in terms of social benefits from the parsimonious use of insufficient public funds in the control of epidemics. Moreover, the Pontryagin's Maximum Principle [17] and the time dependent control strategies are also applied for the studies of HIV models [18] , mosquito-borne diseases [19] , insect transmitted diseases [20] , dengue models [21] , tuberculosis models [22] and others [23] .

In this paper, we consider the optimal and sub-optimal control strategies associated with quarantining asymptomatic individuals and isolating symptomatic individuals for a SARS SEQIJR model presented by Gumel [3] . The model monitors the dynamics of six sub-populations (classes), namely susceptible, asymptomatic, quarantined, symptomatic, isolated, and recovered individuals. We introduce, for the first time, into this model two control variables representing the rate of quarantining of asymptomatic individuals who have been exposed to the virus, but have not yet developed clinical symptoms, and the rate of isolating of symptomatic individuals, respectively.

The rest of the paper is organized as follows. Section 2 describes the SEQIJR SARS model. The analysis of optimization problems and the sub-optimal solution are presented in Sections 3 and 4, respectively. In Section 5, we present the numerical method and the simulation results. Finally, the conclusions are summarized in Section 6.

In this section, we use the dynamic SEQIJR model for SARS presented by Gumel [3] . By choosing the proper parameters, the simulation results from such a model fit the available data fairly well (see Section 3 in [3] ). 

The model synthesizes some demographic effects by assuming a proportional natural death rate µ > 0 in each of the six sub-populations of the model and the constant recruitment rate Λ. In addition, it also includes a net inflow of asymptomatic individuals into the region at a rate p per unit time. This parameter includes new births, immigration and emigration. Owing to the continuing inflow of travellers into the community during the SARS outbreaks, and the absence of a rapid and effective screening test, it is assumed that some of these travellers are asymptomatically infected and enter the community at a rate p per day. In this work, we set p to zero for simplicity. We assume that a susceptible individual may be infected through contacts with an asymptomatic individual, a quarantined individual, an infected individual, or an isolated individual.

The transmission coefficients for these four classes of infected individuals are β, ε E β, ε Q β, and ε J β respectively. Although SARS is believed to be transmitted exclusively by symptomatic individuals, a very low rate of transmission by asymptomatic individuals cannot yet be ruled out. In the model this possibility is depicted by employing the modification parameter ε E , where 0 ε E < 1.

The modification parameter ε Q 0 represents varying levels of hygiene precautions during quarantine, and a similar interpretation is given to the modification parameter ε J during isolation [3] . An asymptomatic individual is transferred into the symptomatic class (at a rate k 1 ), and a quarantined individual is transferred into the isolated class (at a rate k 2 ) because of the development of clinical symptoms. d 1 and d 2 are per-capita disease induced death rates for the symptomatic individuals and isolated individuals, respectively. In addition, σ 1 and σ 2 are per-capita recovery rates for the symptomatic individuals and isolated individuals, respectively.

The control variables, u 1 (t) and u 2 (t), are bounded, Lebesgue integrable functions [22] . The control u 1 (t) represents the rate of quarantining of people who have been in contact with an infected individual by a quarantine program and educational campaigns. The control u 2 (t) represents the rate of isolating of symptomatic individuals by an isolation program. Furthermore, from the epidemiological modeling viewpoint, the transfer rate u 1 E corresponds to an exponential waiting time e −u 1 t , as the fraction that is still in the asymptomatic class t units after entering this class and to 1/u 1 as the mean waiting time [24] . And the interpretation of u 2 I is similar to u 1 E.

The problem is to minimize the cost function

subject to the differential equations (1), where t f is the final time. This performance specification involves the numbers of individuals of symptomatic, asymptomatic, quarantined, or isolated, respectively, as well as the cost for applying quarantine control (u 1 ) and isolation control (u 2 ). The total cost includes not only the consumption for every individual, but also the cost of organization, management, and cooperation etc. Hence, the cost function should be nonlinear. In this paper, a quadratic function is implemented for measuring the control cost by reference to many papers in epidemic control [18] [19] [20] [21] [22] . The coefficients,B 1 , B 2 , B 3 , B 4 , C 1 and C 2 , are balancing cost factors due to size and importance of the six parts of the objective function. We seek to find an optimal control pair, u * 1 and u * 2 , such that

Here

Pontryagin's Maximum Principle [17] provides the necessary conditions for an optimal control problem. This principle converts (1)-(3) into a problem of minimizing a Hamiltonian, H , pointwise with respect to u 1 and u 2 :

where f i is the right hand side of the differential equation of i-th state variable. By applying Pontryagin's Maximum Principle [17] , we have the adjoint equations

with transversality conditions

By the bounds in Ω = {(u 1 , u 2 ) ∈ L 1 (0, t f )|a i u i b i , i = 1, 2}, the optimal control is given by

and

which are derived from the condition

Remark 1. Due to the a priori boundedness of the state and adjoint functions and the resulting Lipschitz structure of the ODEs, we obtain the uniqueness of the optimal control for small t f [22] . The uniqueness of the optimal control pair follows from the uniqueness of the optimality system, which consists of (1), (5) and (6) with characterizations (7) and (8) . There is a restriction on the length of the time interval in order to guarantee the uniqueness of the optimality system. This smallness restriction on the length on the time interval is due to the opposite time orientations of (1), (5) and (6); the state problem has initial values and the adjoint problem has final values. This restriction is very common in control problems (see [18, 22] ).

Remark 2. The problem described above is a Two Point Boundary Value Problem (TPBVP), with specified initial conditions for the state equation (1), and terminal boundary conditions (6) for the adjoint equation (5) . It can be numerically solved by the iteration method presented in Section 5.

In this section, the sub-optimal solution is obtained by a genetic algorithm which follows the idea of [21] . Let the class of admissible controls be restricted to the collection of functions of type

where u i k ∈ Ω are constants, and I i are closed intervals [t i , t i + ∆ i ] such that I i ∩ I j = ∅ if i = j. Therefore, the restricted class admissible controls consist of pulses of height u i k (k = 1, 2) and width ∆ i , starting at time t i . In order to simplify the actual policy implementation, we assume that the two controls are to be of pulsed form. Let u k (t) = u N k (t) (k = 1, 2) be admissible controls with N pulses, characterized by

. Therefore, the sub-optimal control problem is to minimize

where J [·, ·] is the same as defined in (2).

for the case of N pulses, and the sub-optimal control problem is simply

In this work, the search for a global minimum of J (π ) is implemented by an extended method based on the Genetic Algorithm [25] . This extension admits the constraints t i + ∆ i = t i+1 , t N + ∆ N = t f , t i 0, ∆ i 0 and a k u i k b k ,k = 1, 2.

In this section, we investigate the numerical algorithm for optimal control in SARS model. The optimal Quarantine and Isolation Control is obtained by solving the optimality conditions, consisting of twelve ODEs from the state and their adjoint equations. This problem (TPBVP) is solved by using the multiple shooting method (see, for instance, [21, 26, 27] ) as follows:

(1) Select n nodes to divide the interval [0, t f ], or, 0 = t 1 t 2 · · · t n = t f .

(2) Choose y i := [x(t i ), λ(t i )], i = 1, . . . , n.

(3) Integrate (1) and (5) for each time interval [t i , t i+1 ) using y i as the initial conditions, and obtaining y(

If |h| is sufficiently close to a minimum value, then stop; Else modify the initial conditions y i for the next iteration (using, for instance, the Newton-Raphson algorithm to make h i := y new i − y(t i ) = 0) and go to step (3). In this work, the optimal control was computed using the multiple shooting method with 30 nodes, up to a tolerance of 10 −5 for |h|. Considering the four weight factors associated with E, Q, I and J , we choose B 1 = B 2 = B 3 = B 4 = 1, C 1 = 300 and C 2 = 600 to illustrate the optimal strategies. Remark 3. We assume that the weight factor C 2 associated with control u 2 is much larger than C 1 , which is associated with control u 1 . This assumption is based on the following facts: The cost associated with u 1 mainly includes the cost of monitoring and quarantining programs, while the cost associated with u 2 includes the cost of monitoring and isolating programs, and the hospital treatment resource. The ideal weights are very hard to obtain in practice. It needs a lot of work on data mining, analyzing, and fitting. Hence, the acquisition of appropriate practical weights is a difficult problem and it still remains for further investigations. It should be pointed out that the weights in the simulations here are of only theoretical interest to illustrate the control strategies proposed in this paper. 1 . The optimal and sub-optimal control laws.

In the simulation, we study the impact of optimal quarantine and isolation control on SARS during the outbreaks in Beijing. We choose the upper bound of u 1 equals 0.50, according to the reasonable case in China that it took at least an average of 2 days to quarantine the asymptomatic individuals who had been exposed to the virus, but had not yet developed clinical symptoms. And the choice of upper bound of u 2 is similar to u 1 . The sub-optimal control was computed for the case N = 3 and specified time intervals: ∆ 1 = 60, ∆ 2 = 240 and ∆ 3 = 60. Parameters such as population size for the genetic algorithm, generation gap and others were adjusted ad hoc [21] . The data and parameters we used, as shown in Table 1 , are mainly from the published data [3, 9] . In order to evaluate the effect of the uncertainty on the model parameters, simulations were repeated with values of β, ε E , ε Q , ε J , µ, k 1 , k 2 , d 1 , d 2 , σ 1 and σ 2 , altered by random values in the range ±10% (normal distribution), while the inputs u 1 and u 2 were kept identical to the computed optimal values. Fig. 1 shows the time dependent optimal and sub-optimal quarantine and isolation control laws. In order to minimize the total infected individuals, E + Q + I + J , the optimal control u 1 is at its upper bound during more than the half of the year, and then u 1 is steadily decreased to the lower bound, while the optimal control u 2 stays at its upper bound for a short time, about 70 days and then steadily decreases to the lower bound over the rest simulated time. The sub-optimal control laws u 1 and u 2 stay at the upper bound during the first two time intervals and the first time interval, respectively.

In fact, at the beginning of simulated time, both the optimal control and sub-optimal control are staying at their upper bounds in order to quarantine and isolate as many as asymptomatic individuals (E) and symptomatic individuals (I ) to prevent the increasing of the number of the infected individuals. For both u 1 and u 2 , the steady decrease to the lower bound in the optimal control case, and the staying at a proper value in the sub-optimal case, are determined by the balance between the cost of the infected individuals and the cost of the controls. A comparison of the optimal control, the sub-optimal control and the constant controls u 1 ≡ 0.2 and u 2 ≡ 0.2 throughout the simulated time are also implemented (Fig. 2) . It is easy to see that the optimal control and sub-optimal control are much more effective for reducing the number of infected individuals and decreasing the time-span of the epidemic. As normally expected, in the early phase of the epidemic breakouts, keeping the quarantine and isolation controls at their upper bounds will directly lead to decreasing of the number of the infected people. Moreover, Fig. 3 shows that if we take the lower bound values for each control throughout the simulated time, the number of infected individuals would reach about 4.1919 × 10 6 (nearly 35% of the total population). This illustrates that the quarantine and isolation strategies are critically important to control the outbreaks of SARS. In addition, in order to illustrate the overall of picture of the epidemic, the number of infected and recovered individuals under the optimal control, the sub-optimal control, the constant control and the lower bound control are shown, respectively, in Figs. 4-9. The costs achieved by the optimal control and sub-optimal control are much less than the constant control and the lower bound control, as shown in Figs. 10 and 11. It is interesting to note that the cost achieved by the sub-optimal control is very close to that of the optimal control, despite the fact that the sub-optimal control is much easier to implement than the optimal control.

In this paper, we introduced, for the first time, two control variables to represent the quarantine and isolation strategies in epidemic control. The optimal and sub-optimal quarantine and isolation controls for a SEQIJR SARS model were studied. The numerical results show both the optimal and sub-optimal strategy for SARS epidemic control are effective. Furthermore, this idea of introducing these two new control variables is helpful for further discussions on catastrophic epidemic control. It should be pointed out that, as was mentioned in [23] , the ideal time-varying optimal strategy might not be applied in practice easily. Nevertheless, it does provide a reference basis on which to design the practical quasi-optimal control strategies or policies, and to assess their effectiveness. However, the proposed sub- Fig. 10 . Costs of the controls and the disease under the optimal control, sub-optimal control and constant control. optimal control in this paper can provide much simpler strategies, which, as verified by numerical simulations, can also lead to a performance very similar to that achieved by the optimal control. In practical implementations, it is considerably easier to use such a sub-optimal control.

Wresting SARS from uncertainty

Consensus document on the epidemiology of severe acute respiratory syndrome (SARS)

Modelling strategies for controlling SARS outbreaks

SARS outbreaks in Ontario, Hong Kong and Singapore: The role of diagnosis and isolation as a control mechanism

Modeling the SARS epidemic

Transmission dynamics of etiological agent of SARS in Hong Kong: Impact of public health interventions

Transmission dynamics and control of severe acute respiratory syndrome

Epidemiological determinants of spread of causal agent of SARS in Hong Kong

Simulating the SARS outbreak in Beijing with limited data

Model parameters and outbreak control for SARS

Stochastic dynamic model of SARS spreading

Curtailing transmission of severe acute respiratory syndrome within a community and its hospital

Transmission of severe acute respiratory syndrome in dynamical small-world networks

Small world and scale free model of transmission of SARS

A discrete epidemic model for SARS transmission and control in China

Applications of Optimal Control Theory in Biomedicine, Marcel Dekker

The Mathematical Theory of Optimal Processes

Optimal control of the chemotherapy of HIV

Optimization of investment policies in the control of mosquito-borne diseases

Optimization of the investment in educational campaigns in the control of insect transmited diseases

Optimal and sub-optimal control in Dengue epidemics

Optimal control of treatments in a two-strain tuberculosis model

The mathematics of infectious diseases

Genetic Algorithms and Simulated Annealing

Introduction to Numerical Analysis

Real time computation of feedback controls for constrained optimal control problems. Part 2: A correction method based on multiple shooting

This work was jointly supported by the National Natural Science Foundation of China under Grant Nos. 60474078, 60574015, 60304001.