key: cord-0542093-x1a6xs5o
authors: Hamelin, Fr'ed'eric; Iggidr, Abderrahman; Rapaport, Alain; Sallet, Gauthier; Souza, Max
title: Identifiability and observability of the SIR model with quarantine
date: 2021-03-08
journal: nan
DOI: nan
sha: 68efc5df3751432fbb7d2ff9eb5a1dc51fd3cc8b
doc_id: 542093
cord_uid: x1a6xs5o

We analyze the identifiability and observability of the well-known SIR epidemic model with an additional compartment Q of the sub-population of infected individuals that are placed in quarantine (SIQR model), considering that the flow of individuals placed in quarantine and the size of the quarantine population are known at any time. Then, we focus on the problem of identification of the model parameters, with the synthesis of an observer.

Many papers in epidemiology proposing a mathematical model using dynamical systems, face the problem of parameter estimation. In general, some parameters are given, extracted from the literature, while the remaining unknown parameters are estimated by fitting the model to some observed data, usually by means of an optimization algorithm based on least squares or maximum likelihood methods. Nevertheless, relatively few studies about the intrinsic property of a model to admit a unique set of parameters values for a given choice of measured variables. On the other hand, this a question that is well known in automatic control. Investigation of identifiability in mathematical epidemiology is relatively recent [11, 8, 13, 10, 6, 1, 2] . Indeed, to the best of our knowledge, the first paper considering the problem in a epidemic model studies identifiability of an intra-host model of HIV, and it has been published in 2003 in an Automatic Control journal [13] . It is also surprising that the observability and identifiability of the original Kermack-Mckendrick model has not been more studied, since this system has been widely used to model an outbreak of an infection. The observability and identifiability of the classical model SIR, with demography and constant population, has been first studied in 2005 [2] .

The global health crisis of COVID-19 outbreak has led to a spectacular resurgence of interest in this type of models, but with specificities related to the detection and isolation of infected individuals [9] . This is why we revisit here the issues of identification and observability for an extended 'SIQR' model [4] for which such an analysis had not yet been performed (to the best of our knowledge). Once the question of identification has been settled, we also tackle the task of proposing a practical strategy for reconstructing the unique values of the parameters.

Inspired by [4, 7] , we consider the classical SIR model (see for instance [5] ), where S, I, R denote the size of the populations of respectively susceptible, infected and recovered individuals, with an additive compartment where Q denotes the size of the population of identified and isolated infectious individuals that have been removed from the infected population and placed in quarantine:

When the size of the total population N is large and the size of the population placed in quarantine remains small compared to N during the considered interval of time, one can consider a simplified model:

(2)

Note that for both models, one has

These models have three parameters: the infectivity parameter β , the recovery rate ρ, that we assume to be identical for the infected populations placed in quarantine or not, and the rate of placement in quarantine α. Theses parameters are unknown but we assume the following hypothesis. 

This assumption implies that the epidemic can spread in the population i.e. at initial time with S(0) = N − I(0) close to N one hasİ(0) > 0.

We assume that • the flow αI(t) of infected people placed in quarantine is known at any time t ≥ 0

• the size Q(t) of the population placed in quarantine is perfectly known at any time t ≥ 0

• the size N of the total population is known • at initial time 0, one has S(0) = N − ε, I(0) = ε, Q(0) = 0, R(0) = 0 with ε ∈ (0, N).

We consider then the observation function

and follow the usual definitions of identifiability and observability of systems [12, 3] . However, note that when Q = 0, the system is not infinitesimally identifiable: the knowledge of the outputs and all its derivative do not allow to determine formally ρ. At I = 0, the system is not identifiable neither. We adopt the following definition of identifiability for these models.

Definition 1. Given N > 0 and ε ∈ (0, N), we shall say that system

is injective, where (S(·), I(·), Q(·), R(·)) is solution of the Cauchy problem for the differential system (1) resp. (2) with

is injective, then the system (1) resp. (2) is identifiable and observable for the observation (3).

4 Analysis of the first model Proof. It consists in showing that parameters and unmeasured variables S and I can be expressed as functions of the successive derivatives of the output vector y. As the variable I cannot reach 0 in finite time, we shall assume I = 0 in the following.

Note first that with Q(0) = 0 one hasQ(0) > 0 and then y 2 (t) = Q(t) > 0 for any t > 0. The dynamics of Q gives directly the expression of the parameter ρ as:

Posit h 1 :=ẏ 1 y 1 . One has from the dynamics of I

and then

Using the equality β S N − Q = α + h 1 + ρ from (5), one obtains

Let us write the derivative of h 2 :ḣ 2 =ḣ 1 (−β I +Q)

which can be also expressed asḣ

Then, using relation (7), one obtains the expressioṅ

or equivalently the equationḣ

to be fulfilled.

Observe that this last equation is a second order polynomial in the variable X = −β I +Q. From (6) and R 0 > 1 one has

and this allows us to show that one also has

Observe also that one has X(0) = −β I(0) +Q(t) = ε (−β + α) < 0. Therefore, by continuity w.r.t. t, we obtain that for t > 0 small enough, X is the unique negative solution oḟ

that is

The parameter α can be then obtained from equation (6) α

where ρ is given by (4) . The initial condition ε is simply reconstructed by ε = y 1 (0)/α and finally one obtains the parameter β = α − X(0)/ε. Proof. As for model (1), one can determine the parameter ρ from any positive time as

Then from the dynamics of I one can write

where h 1 is a known function. Differentiating h 1 with respect to the time giveṡ

and differentiating twiceḧ

With the expression (10), we rewrite this last equation as follows

and with expression (11)ḧ

One obtains then the following expression for β I(t)

(note that this expression is well defined because h 1 (t) + α > 0 for any t and thusḣ 1 (t) < 0).

Finally, from (11) one reconstructs the parameter

and then the parameter

At last, the initial condition is recovered as ε = y 1 (0)/α.

The former analyses have shown that models (1) and (2) are not infinitesimally identifiable at time 0 when the initial state is (N − ε, ε, 0, 0). One has to wait a short time t > 0 to have Q(t) > 0 and formally identify parameters. Thus, we expect a weak accuracy of the parameters estimation at the very beginning, that should improve with time while the state get away from this initial state and new measurements come. This is why we have opted for a dynamical estimation with the help of observers. The use of observers, although usually dedicated to state estimation (rather to parameters estimation) possesses the advantage to tune the speed of error decay. Moreover, the choice of a speed-accuracy compromise can be balanced thru simulations with synthetic data corrupted by noise, before facing real data. Note also that for large times, the solutions of (1) and (2) converge asymptotically to non-observable non-identifiable states when I and Q are both null. Consequently, we do not look precisely for results about asymptotic convergence of the error (as this is usually done in observers theory), but rather for an exponential decay of the error estimation during initial transients.

In this section, we shall consider the additional hypothesis Assumption 2. One has α ≤ ρ that means that the rate of placement in quarantine is not larger than the recovery rate, which is often the case in epidemic regimes.

We shall denote the elementary symmetric polynomials, where X is a vector in R n , as

Proposition 3. Let λ and µ be two positive vectors in R 4 and R 3 respectively. For t > 0, consider the dynamical system

) − (ẑ 2 − log(y 2 (t))),

with the gains vector

is an estimator of [ρ, β , α] , whose exponential decay of the error can be made as fast as desired keeping the number

large, as long as S(t)/N remains close to 1. Moreover, the state I is estimated witĥ

Proof. Posit δ := β − α.

As far as S/N remains close to 1, the size of the population I is small compared to S, and the dynamics of the outputs y 1 = αI and y 2 = Q can be approximated by the linear dynamics ẏ 1 = δ y 1 − ρy 1

For t > 0 and i = 1, 2, consider the new outputs z i (t) = log(y i (t)), whose dynamics is given by the system

For this sub-system with unknown parameters δ and ρ, we consider the following candidate observer in R 4 :

The dynamics of the error e 1 = (ẑ 1 ,ẑ 2 ,δ ,ρ) − (z 1 , z 2 , δ , ρ) is given by the linear time invariant systemė 1 = M 1 e 1 with

A calculation shows that with the choice

. This shows that the first four equations of system (12) gives the reconstruction of the parameters δ and ρ with an exponential decay of the error larger than min i λ i .

and consider the variable

As far as S(T )/N remains close to 1, one can write the approximatioṅ

Then, this amounts to approximate the dynamics of (1) or (2) in the (y 1 , v, k) coordinates by the following dynamical system

where k is an unknown parameter, for which we consider the following candidate observer in R 3

whose dynamics of the error e 2 = (ŷ 1 ,v,k) − (y 1 , v, k) is given by the non-autonomous linear systeṁ

with

One can easily check that for the choice

Then, from (22), we obtain the upper bound on the error decrease

T 0 y 1 (τ)dτ ||e 2 (0)|| whose exponential decay can be made as large as desired with large l.

Finally, from the reconstruction of parameters δ , k by observers (18), (21) and expressions (15) and (19), the original parameters α, β are recovered as roots of

Note first that Assumption 1 implies β > α and thus k 2 − 4kδ > 0. Moreover, one has k > β (1 + ρ α ), and by Assumption 2 one has k > 2β , which implies that only the smaller root of (23) is valid, leading to the expression (13) of the estimator. Note that this expression preserves the exponential decay of the error obtained for δ , k and ρ.

Let us make some comments about this observer. It consists in reconstructing functions of the parameters δ and k and not directly the parameters α, β . There is an apparent redundancy of variablesẑ 1 andŷ 1 in dynamics (12) , which reconstruct log y 1 and y 1 . Indeed, this allows to decouple the observer into two sub-systems of dimensions 4 and 3, which avoids the use of two large correction gains compared to a full order observer. Finally, outputs of these two sub-systems are coupled in expression (13) to reconstruct the original parameters.

The proposed observer has been tested with synthetic data for a population size N = 10 5 with parameter values α = 0.07, β = 0.4, ρ = 0.1 and initial condition I(0) = 10, Q(0) = 5, R(0) = 0 over a time horizon of 10 days (see Figure 1 ). The gains have been computed for the choice of vectors λ = [1; 1.5; 2; 2.5] and µ = [1/(13.10 3 ); 1/(15.10 3 ); 1/(19.10 3 )]. Note that vector µ has been chosen quite small to avoid too large gains when multiplied by N in the observer equations.

Then, we have simulated a measurement noise with a centered Gaussian law of variance equal to 5% of the signal (see Figure 2 ). Because the expression (14) of the estimation of the state I is not filtered, we have applied a moving average smoothing to the estimationÎ (see Figure 3) .

Finally, these simulations show that the method allows to reconstruct the parameter values in few days in a quite accurate manner. The estimation of the size of the infected population I allows then the use of the model for predictions of the epidemics.

This work shows that although the identifiability of the SIR-Q models has singularity points where measured variables are null, it is possible to design an observer with exponential decay of the estimation error during the first stage of the epidemics, and recover parameters in few days. Further investigations will concern real data of COVID epidemics provided by various territories. 

Identifiability and estimation of multiple transmission pathways in Cholera and waterborne disease

The structural identifiability of the susceptible infected recovered model with seasonal forcing

Observability for systems with more outputs than inputs and asymptotic observers

Effects of quarantine in six endemic models for infectious diseases

An introduction to mathematical modeling of infectious diseases

On identifiability of nonlinear ODE models and applications in viral dynamics

Mathematical models of influenza: the role of crossimmunity, quarantine and age-structure

Identifiability analysis of an epidemiological model in a structured population

Why is it difficult to accurately predict the covid-19 epidemic?

An effective automatic procedure for testing parameter identifiability of HIV/AIDS models

Structural and practical identifiability issues of immunoepidemiological vector-host models with application to Rift Valley Fever

From experimental data, Translated from the 1994 French original and revised by the authors

Identifiability of nonlinear systems with application to HIV/AIDS models