key: cord-0679594-2sjfhemf
authors: Batistela, Cristiane M.; Correa, Diego P. F.; Bueno, 'Atila M; Piqueira, Jos'e R. C.
title: Compartmental model with loss of immunity: analysis and parameters estimation for Covid-19
date: 2020-07-02
journal: nan
DOI: nan
sha: 23fd470a46179ae76ee55a3818c68ad20fdcf06b
doc_id: 679594
cord_uid: 2sjfhemf

The outbreak of Covid-19 led the world to an unprecedent health and economical crisis. In an attempt to responde to this emergency researchers worldwide are intensively studying the Covid-19 pandemic dynamics. In this work, a SIRSi compartmental model is proposed, which is a modification of the known classical SIR model. The proposed SIRSi model considers differences in the immunization within a population, and the possibility of unreported or asymptomatic cases. The model is adjusted to three major cities of S~ao Paulo State, in Brazil, namely, S~ao Paulo, Santos and Campinas, providing estimates on the duration and peaks of the outbreak.

The Wuhan Municipal Health Commission reported a cluster of 27 pneumonia cases on 31st December 2019, in the Wuhan city, Hubei Province in China.

On 1st January 2020 the World Health Organization (WHO) set up the Incident Management Support Team putting the organization to an emergency level for dealing with the outbreak. On 5th January 2020 WHO published the first outbreak news on the new virus. On 7th January the causative agent was identified and named Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2) (WHO named the disease Covid-19). On 13th January 2020 the first case outside China was reported in Thailand. On 22nd January 2020 WHO stated that there was evidence of human-to-human transmission, and approximately seven weeks later, on 13th March 2020, WHO characterized the Covid-19 as a pandemic [1, 2] .

In Brazil, the first confirmed Covid-19 case was reported on the 26th February 2020, and up to 25th June 2020 there was 1,188,631 confirmed cases with 53,830 deaths [3] . Globally, according to [4] there was 9,292,202 confirmed cases, and 479,133 deaths up to 25th June 2020 1 .

Most cases are asymptomatic carriers and are spontaneously resolved, however, some developed various fatal complications, notably for patients with comorbidities [2] , and, allied to the fastly spread of covid-19, the worldwide emergency state brings up with important, and yet unanswered, questions related to the contagion dynamical behavior and its mitigation and control strategies. As a result, strategies to contain the contagion such as social distancing, quarantine and complete lockdown of areas have been studied [5, 6, 7, 8, 9, 10] In Mexico, on 18th March 2020, the Mexican Health Secretariat reported that the pandemic was going to last 90 days, with 250,656 expected cases. On the next day, the Health Secretariat informed that approximately 9.8% (24, 594) of the cases would need hospitalization, and 4.2% (10,528) of the total cases would be critical patients, needing intensive care. On the same date the number of available Health Care units at that time was 4,291 with 2,053 ventilators [11] .

The situation in Mexico led to the adoption of non-pharmaceutical interventions, such as washing hands, social distancing, cough/sneeze etiquette, and so on.

In Rio de Janeiro, Brazil, the social distancing started on 17 March 2020 [8] , and the government of the state of São Paulo, Brazil, decreed quarantine on 1 Reported on 2:25pm CEST.

23rd March 2020. Up to 25th June 2020, the State of São Paulo had 248,587 confirmed cases with 13,759 deaths.

One important question is related to the patient's immunity after recovery.

The difference in immunity after recovery have been reported in humans in [12] , and in experiments with rhesus macaques [13] . The experiment in [12] collected plasma from 175 Covid-19 recovered patients, and SARS-CoV-2-specific neutralizing antibodies (NAbs) were detected from day 10 to 15 after the onset of the disease and remained thereafter. Nonetheless, the NAbs levels were variable in the cohort, 52 (≈30%) of the patients developed low levels of these antibodies, from which with 10 (≈6% of the recovered patients) the NAbs levels were undetectable, 25 (≈14%) developed high levels.

In [14] the data suggests that the response to SARS-CoV-2 is imbalanced regarding to controlling the virus replication versus the activation of the adaptive immune response. In some cases the immune system doesn't work properly and lung cells remain vulnerable to infection. The virus continues replicating while the immune response system attacks infected cells, killing even healthy nearby cells, and the lung tissue becomes seriously inflamed. This seems to be the mechanism that make some patients become severely ill weeks after their initial infection. Additionally, SARS-CoV-2 probably induces immunity like other coronaviruses, however, it is so new that this mechanism isn't fully understood [15] .

The economic crisis due to the pandemics is another important issue. According to [16] the mortality rate of Covid-19 is not necessarily correlated with the economic risk to global economy since governments, companies, consumers and media reacted to the economic shock. However, a global recession seems to be inevitable, its duration and intensity will depende on the success of measures to prevent the spread of Covid-19.

Taking the whole scenario into account researchers worldwide are intensively studying and developing mathematical models of the Covid-19 outbreak. The knowledge on this pandemic dynamics is important for providing estimates on the duration and peaks of the outbreak.

The macro-modelling of infectious disease spread have been a field of research from many years since the simple deterministic model of Kermack and McKendrick [17, 18, 19, 20] , that provides a dynamical model classifying individuals in a population as Susceptible -Infected -Removed (SIR) [21] . A number of classical models such as SIR [22] and Susceptible -Exposed -Infected -Removed (SEIR) [23, 24] have been proposed for epidemic modelling.

In addition, works considering time delayed [25] and fractional order [26] dynamical systems applied to Covid-19 outbreak have also been proposed.

In [27] a SIRD model had been adjusted to Covid-19 spreads in China, Italy and France. The results have shown that the recovery rate were similar for the situations of China, Itally and France. In [28] , a mixed analytical-statistical inverse problem is used to predict the Covid-29 progression in Brazil, a SIRU model, where U stands for Unreported cases, was used for the direct problem, and a mixture of parameter estimation with Bayesian inference for the inverse problem analysis.

Compartmental models considering immigration and home isolation, or quarantine, are proposed on [10] , all the situations presented infection-endemic equilibrium, the results demonstrated that home isolation, or lockdown, mitigates the chance of infection.

In this work the proposed model is a modification of the original compartmental SIR model of Kermack and McKendrick [17, 18, 19] , including a sick (S ick ) population compartment representing nodes of the network that manifest the symptoms of the disease. The proposed Susceptible -Infected -Removed -Sick (SIRSi) model also considers the birth and death of individuals in the given population. In addition, a feedback from those recovered who did not gain immunity or loss their immunity after a period of time and become susceptible again is also introduced. 

The proposed SIRSi model can be seen in Fig. 1 . In this model, the susceptible population S is infected at a rate α when contact infected individuals from I. The susceptible compartment also receives a population, at a rate 1/γ, who didn't gain complete immunity after recovering or who loss their immunity after a period of time.

The compartment I represents the infectious population in incubation stage prior the onset of symptoms. Infection transmission during this period has been reported in [29, 30, 31, 32] . The infected population can be asymptomatic or symptomatic. The period between infection and onset of symptoms, 1/β 2 , ranges from 3 to 38 days with median of 5.2. Once the infected individual is tested positive and the case is documented, the case is moved to the S ick compartment. Those who don't develop severe symptoms become asymptomatic.

In [33] the estimative of the infections originated from undocumented cases is as high as 86%, which include mild, limited and lack of symptoms infectious individuals. In other recent studies [34, 35] , it was found that, 20% -40% of positive tested patients were asymptomatic. In this work, and according to [36] , the asymptomatic population is considered as those of under-reporting cases. This population could be as 7 times bigger than the size of the documented cases. This group, under-reporting cases, become recovered with the period 1/β 1 .

Some of the individuals within this population could eventually develop symptoms. In [37] , it has been reported that the average time between infection and the onset of symptoms can be 4.6 days. Once the case becomes documented it should be moved to the S ick compartment. S ick are those infected which seek for medical attention with severe symptoms. This population are those who tested positive for COVID-19. In [38] it was reported that this population could represent up to 19.9% of the total documented cases, of which 13.8% are severe cases and 6.1% require intensive care. The sick population become recovered with the period 1/β 3 or removed at a rate σ (see Fig. 1 ).

In order to consider the social distancing measure effect in the number of infected and deaths toll shown in Fig. 1 , the parameter θ is introduced on the mathematical model (1), where θ is subject to the constraint 0 < θ < 1.

Consequently, given these facts, the mathematical SIRSi model is given by (1),

where λ and δ are the birth and death rates, respectively.

It is important to notice that the number of documented cases is a key information that should emerge, somehow, from (1). The reason is the fact that the accumulated number of confirmed cases is publicly available, and will be used to fine-tuning the model. The number of daily new infections is also available, but it tends to be noisy and less representative of the dynamics. (1), and parameters α, β 1 , β 2 , β 3 , σ, γ, θ ∈ R + .

To investigate the influence of the introduction of the feedback from those recovered with no immunity which become susceptible again, and dividing the population into groups, the equilibrium points related to the model described by (1) must be determined and their stability discussed.

Assuming α = 0, i.e., susceptible can be converted into infected, and despite being an assumption it is realistic for a spreading disease, the equilibrium points are such that f (x * ) = 0.

Using the Hartman-Grobman Theorem [39] the local stability of the equilibrium points can be determined by the eigenvalues of the Jacobian matrix computed on each equilibrium point. The Jacobian J = ∂f ∂x x * of (1) is given by (2) .

In the following the disease-free (section 3.1) and the endemic (section 3.2) equilibrium points are determined.

The disease-free equilibrium is a state corresponding to the absence of infected individuals, i.e., I * = 0. Applying this condition to the equilibrium of (1), the point can be determined.

Assume that there exist a disease free-equilibrium (

Considering P 1 , the corresponding linear system Jacobian [39] J P1 = ∂f ∂x x * is given by (4) .

By the Laplace determinant development [40] , the eigenvalues of (4) are the elements in the diagonal, that is,

The eigenvalues above with the Hartman-Grobman Theorem indicates that

(1) presents one disease-free equilibrium point, and considering the condition given by the eigenvalue

the eigenvector associated indicates an asymptotically stable direction. Con-

indicating a bifurcation in the parameter space.

The endemic equilibrium points are characterized by the existence of infected people in the population, that is, (I * = 0).

Therefore, assuming the existence of an endemic equilibrium point, with

x * ∈ U, such that f (x * ) = 0, the equilibrium point P 2 = x * = (S * , I * , R * , S * ick ) T in the first octant of R 4 is given by (6),

Accordingly, the existence condition of a positive endemic equilibrium P 2 =

x * = (S * , I * , R * , S * ick ) T ∈ R + is given by (6) .

Note that condition (6) reflects the fact that, in order to the endemic equilibrium exists, the rate α at which, the people flux rate, represented by λ/δ, become infected, has to be greater than the rate at which infected population leave the compartment I, either overcoming the disease or becoming symptomatic.

It is important to highlight that λ/δ could also represent the total people flux commuting from a different city in a multi-population model.

The linearization A = J P2 = ∂f ∂x x * at the endemic equilibrium is:

The characteristic polynomial det(A − I d ξ) = 0 is

with

Any further effort to analytically analyze ξ eigenvalues, becomes quite difficult due to the coefficients complexity. A possible alternative approach is to go for numerical calculations.

However, some insight for the model with feedback γ can be obtained, in terms of bifurcations and stability, if we analyze the eigenvalues when γ = 0 and γ = 0.

Note that in this case, the endemic equilibrium is still possible. Computing the eigenvalues results,

The eigenvalues ξ 3 and ξ 4 can be either complex conjugate stable, or both real. The eigenvalue ξ 3 needs to be further studied due to the possibility of bifurcation.

Analysing the eigenvalue ξ 3 , if α(1 − θ) λ > √ ∆, P 2 is stable, and consequently (6) holds true.

On the other hand, if α(1 − θ) λ < √ ∆, the endemic equilibrium point is unstable, and the existence condition (6) is not satisfied, consequently, the endemic equilibrium point P 2 , if existing, is stable.

Another insight can be obtained by looking with the case γ → ∞. In this case, the endemic equilibrium becomes:

which is subject to the same existence condition given in equation (6) .

Under the assumption γ → ∞, the characteristic polynomials coefficients in equation (9) can be approximated by (12):

the characteristic polynomial (8) can be rewritten as in (13):

assuming that at least one root |ξ| goes to infinity as γ → ∞, we rewrite the polynomial as

then, ξ 1 = −γ, and γ → ∞, so, one eigenvalue seem to be going to −∞ as γ → ∞. To find an approximation to the other three roots, the characteristic polynomial is rewritten in equation (13) . It can be also assumed that the other three roots are finite,

Looking for insight numerical experiments are performed.

In this section the proposed SIRSi model (1) (see Fig. 1 The data available for the birth and death rates from the public repository was out of date, with the last report been published in 2018, therefore, for 2019 and 2020 these rates were computed using interpolation, results are shown in For the parameters fitting, the feedback parameter γ was kept constant at values γ ∈ {0, 0.01, 0.02, 0.03, 0.04}, and the trust-region reflective least-squares algorithm [43, 44, 45] was used for each one of the three cities into consideration.

All parameters and initial conditions computed are normalized with respect to the total population in each case. The results are shown in tables 2 to 7.

In order to assess the influence of the parameter γ in the endemic equilibrium, the eigenvalues were plotted with the set of fitted parameters for γ from 0 to 2.

If figures 4 are shown the eigenvalues for the endemic equilibrium for each city computed for each one of the fitted sets. At γ = 0 eigenvalues are stable for 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 S ick0 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 R0 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 1.008e+00 1.008e+00 1.007e+00 1.007e+00 1.005e+00 I0 8.749e-06 9.148e-06 9.364e-06 9.614e-06 1.776e-05 S ick0 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 R0 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 In 5 a closer view of the eigenvalues around the origin are shown.

In this section the numerical experiments are conducted using the MATLAB-Simulink [46] for two different conditions. Firstly, the SIRSi model is fitted with the real data for the S icks population, and for different values of the parameter γ. In the sequence, the simulation for the infected population I, that can be inferred from the proposed model, is carried out.

The numerical experiments, as in section 4, were conducted for three major cities in the state of São Paulo, namely, São Paulo, Campinas, and Santos.

The initial condition is x 0 = (S 0 , I 0 , S ick0 , R 0 ) T , where S 0 and I 0 are the normalized susceptible and infected populations, which are considered free parameters in the sense that they can be modified by the fitting algorithm.

In Fig.6 it can be seen that the SIRSi model is adjusted for the confirmed cases of infected people data.

Considering that the acquired immunity is permanent, i.e., γ = 0, and that the isolation rate is constant, the peak of the infection occurs soon after July 2020, and until the end of the the same year, the disease will not persist, since the number of confirmed cases will go down to zero.

On the other hand, assuming that immunity is not permanent and adopting a reinfection rate γ = 0.01, meaning that every 100 days a recovered person becomes susceptible again, the model predicts a decrease in the confirmed cases Confirmed cases in São Paulo city 10 4 Real data =0 =0.01 =0.02 =0.03 =0.04 In Fig.7 , the infected compartment I inferred from the SIRSi model is presented, showing that the peak of infection is close to July 2020.

Increasing γ will reduce the time for a recovered person to become susceptible again, causing the peaks in Fig. 7 to increase, when compared to the curves for lower values of γ. This behavior, however, cannot be observed in Fig.6 , indicating that the increase the re-susceptibility feedback gain γ possibly contributes to the increase of asymptomatic or unreported infected cases.

In addition, it appears that if those recovered acquire permanent immunity γ = 0, the number of infected people tends to zero by the end of 2020. For γ = 0.01, it appears that there is a small increase in January 2022. For γ = 0.02 a new wave of confirmed cases can be seen in Fig.6 , and accordingly, the increase in the infected population is also observed in Fig.7 . For São Paulo, the numerical experiments show that considering any reinfection rate, there will be confirmed infected cases and unreported infected cases until 2023, indicating the need for a control strategy, being necessary and the study of preventive inoculation. people data for Santos. Assuming that the immunity acquired is permanent, γ = 0, and that the isolation rate is constant, the peak of the confirmed cases in Santos will occur very close to July 2020. and similarly to São Paulo (see Fig.   6 ), until the end of the same year, the disease will not persist with the number of confirmed cases going down to zero.

Adopting a nonzero reinfection rate, such that one person every 100 days For γ = 0.04, after the peak of the confirmed cases, a second wave can be observed in the numerical results before the end of 2020, delaying the reduction of the confirmed cases.

For Santos, the numerical experiments show that for γ = 0 and for γ = 0.03, the numbers of confirmed cases tend to zero in the beginning of 2023.

The infected compartment I inferred from the SIRSi model is presented in Fig. 9 , showing that the peak of infection is close to July 2020. The increase in re-susceptibility feedback gain γ, will reduction the time for an infected person to be susceptible again, causing higher peaks when compared with the curves for lower values of γ. This behavior does not occur in Fig.8 indicating that the increase in feedback possibly contributes to the increase in asymptomatic or unreported infected cases. This situation is similar to what is observed for São Paulo.

In addition, for γ = 0, the number of the infected people I tends to zero before the end of the 2020 (see the curve for γ = 0 in Fig.9 ).

For γ = 0.01, a new wave of infection in 2021 can be seen, and for γ = 0.02 the peak of the second wave of infection is near January 2021 (see Fig.9 , γ = 0.01 and γ = 0.02).

Unlike São Paulo, the highest peak of infection among the unreported occurs when γ = 0.03 and this behavior suggests a more detailed study of the dynamics, because together with the situation in which the infected person acquires permanent immunity , these rates suggest that the saving of confirmed cases (see Fig.8 and asymptomatic infected individuals tends to zero more quickly.

In the situation in which a recovered person is liable to a new susceptibility in 25 days, it is observed that the infection persists in the population for a longer time, as shown by the curve with γ = 0.04 in Fig.9 and justifies the policy strategies public policies, including vaccination.

In Fig. 10 the SIRSi model adjusted to the confirmed cases of infected people data in Campinas is shown.

For permanent acquired immunity, γ = 0, and constant isolation rate, the peak of confirmed infection cases occur in the beginning of the second half of 2020.

Considering the re-susceptibility feedback gain γ > 0, in Fig. 10 , it seems that with the increasing of γ the time necessary for the number of confirmed infected cases go down to zero is slightly bigger, unlike the other two cities studied. In addition, Campinas does not present a second wave of infection, even with the variation γ.

The general behavior of Campinas, concerning the sensitivity analysis for γ, present results that differ from Santos and São Paulo. It can be noticed in The infected compartment of Campinas presents the peak of infection close to the beginning of the second half of 2020, Fig.11 .

The increase in the reinfection parameter, causes the peak to increase and this occurs in the figure Fig.10 indicating that the increase in feedback possibly contributes to the increase in asymptomatic or unreported infected.

The proposed SIRSi model was fitted to publicly available data of the Covid-19 outbreak, providing estimates on the duration and peaks of the outbreak.

In addition, the model allows to infer information related to unreported and asymptomatic cases. For the city of Campinas, it is suggested to collect more data, because as the data of the confirmed infected presents a certain distance from the peak of the infection, the dynamics of the model may undergo some significant change,

given the sensitivity of the model to disturbances.

Data are publicly available with [41, 42] .

There is no conflict of interest between the authors.

World Health Organization. WHO Timeline -COVID-19

World health organization declares global emergency: A review of the

Ministério da Saúde do Brasil

World Health Organization. WHO Coronavirus Disease (Covid-19) Dashboard [online

Development of new hybrid model of discrete wavelet decomposition and autoregressive integrated movin Chaos

The effect of control strategies to reduce social mixing on outcomes of the covid-19 epidemic in wuhan

Feasibility of controlling covid-19 outbreaks by isolation of cases and contacts

Covid-19 spreading in rio de janeiro, brazil: Do the policies of social isolation really work?

Strategies to deal with the COVID-19 pandemic, Brazilian Oral Research 34

Covid-19 created chaos across the globe: Three novel quarantine epidemic models

Modeling behavioral change and covid-19 containment in mexico: A trade-off between lockdown and com

Neutralizing antibody responses to sars-cov-2 in a covid-19 recovered patient cohort and their implication

Lack of reinfection in rhesus macaques infected with sars-cov-2, bioRxiv

Imbalanced host response to sars-cov-2 drives development of covid-19

After recovering from covid-19

Economic effects of coronavirus outbreak (covid-19) on the world economy

A contribution to the mathematical theory of epidemics

Contributions to the mathematical theory of epidemics. ii. &#x2014;the problem of endemicity

Contributions to the mathematical theory of epidemics. iii.&#x2014;further studies of the problem of end

Macro-modelling and prediction of epidemic spread at community level

Considering quarantine in the sira malware propagation model

Estimation of parameters in a structured sir model

A double epidemic model for the sars propagation

Seir modeling of the italian epidemic of sars-cov-2 using computational swarm intelligence

On a comprehensive model of the novel coronavirus (covid-19) under mittag-leffler derivative

Analysis and forecast of covid-19 spreading in china, italy and france

Naveira-Cotta, p. magal, Parametric identification and public health measures influence on the covid-19 epidemic evolution in braz medRxiv

Transmission of 2019-NCOV infection from an asymptomatic contact in Germany

Familial cluster of COVID-19 infection from an asymptomatic

A familial cluster of pneumonia associated with the 2019 novel coronavirus indicating person-to-person tr

Characteristics of COVID-19 infection in Beijing

Substantial undocumented infection facilitates the rapid dissemination of novel coronavirus (sars-cov-2)

Estimating the asymptomatic proportion of coronavirus disease 2019 (covid-19) cases on board the diamo

Estimation of the asymptomatic ratio of novel coronavirus infections (covid-19), medRxiv (2020)

Estimate of underreporting of COVID-19 in Brazil by Acute Respiratory S Notas Técnicas Cedeplar-UFMG 010, Cedeplar

Impact of non-pharmaceutical interventions (npis) to reduce covid-19 mortality and healthcare demand

Ferguson, Estimates of the severity of covid-19 disease, medRxiv

Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields

Differential Equations, Dynamical Systems, and Linear Algebra

Portal estatísticas do Estado de São Paulo

On the implementation of a primal-dual interior point method

Recursive quadratic programming methods based on the augmented lagrangian Computation Mathematical Programming

Numerical Computing with MATLAB, Other titles in applied mathematics

A time delay dynamic system with external source for the local outbreak of 2019-ncov,

Analysis 0 (0) (2020) 1-12.This figure "unesp.gif" is available in "gif" format from:http://arxiv.org/ps/2007.01295v2