key: cord-0265490-jgko7b0j
authors: Macedo, A. M. S.; Brum, A. A.; Duarte-Filho, G. C.; Almeida, F. A. G.; Ospina, R.; Vasconcelos, G. L.
title: A comparative analysis between a SIRD compartmental model and the Richards growth model
date: 2020-08-06
journal: nan
DOI: 10.1101/2020.08.04.20168120
sha: e95350fe1b6b44ecf413a8efc67ab0b2fda46ca6
doc_id: 265490
cord_uid: jgko7b0j

We propose a compartmental SIRD model with time-dependent parameters that can be used to give epidemiological interpretations to the phenomenological parameters of the Richards growth model. We illustrate the use of the map between these two models by fitting the fatality curves of the COVID-19 epidemic data in Italy, Germany, Sweden, Netherlands, Cuba, and Japan.

The pandemic of the novel coronavirus disease (COVID-19) has created a major worldwide sanitary crisis [1, 2] . Developing a proper understanding of the dynamics of the COVID-19 epidemic curves is an ongoing challenge. In modeling epidemics, in general, compartmental models [3] have been to some extent the tool of choice. However, in the particular case of the COVID-19 epidemic, standard compartmental models, such as SIR, SEIR, and SIRD, have so far failed to produce a good description of the empirical data, despite a great amount of intensive work [4, 5, 6, 7, 8, 9, 10, 11] . In this context, phenomenological growth models have met with some success, particularly in the description of cumulative death curves [12, 13, 14] . The recent discovery, within the context of a generalized growth model known as the beta logistic model [15] , of a slow, power-law approach towards the plateau in the final stage of the epidemic curves is another remarkable example of this qualitative success. Growth models, however, have the drawback that their parameters may not be easily interpreted in terms of standard epidemiological concepts [16] , as can the parameters of the usual compartmental models.

As a concrete example, consider the transmission rate parameter β of the SIR model [3] . It can be easily interpreted as the probability that a contact between a susceptible individual and an infective one leads to a transmission of the pathogen, times the number of contacts per day. Although the value of β cannot be measured directly in a model independent way, and it is probably not even constant in the COVID-19 epidemic curves, the epidemiological meaning of the parameter is nonetheless easy to grasp conceptually. As a result, models that incorporate such parameters in their basic equations are sometimes regarded as "more epidemiological," so to speak, than others that do not use similar parameters. This state of affairs creates a somewhat paradoxical scenario, in which we have, on the one hand, the striking empirical success of phenomenological growth models sometimes being downplayed, owing to the lack of a simple epidemiological picture of the underlying mechanism [16] , and, on the other hand, the failure of traditional epidemiological compartmental models to produce good quantitative agreement with the empirical COVID-19 data. A glaring instance of the inadequacy of standard compartmental models for the COVID-19 epidemic is their inability to predict the power-law behavior often seen in the early-growth regime as well as in the saturation phase of the accumulated death curves-a feature that is well captured by growth models [15] , as already mentioned.

It is clear that a kind of compromise is highly desirable, in which we get the benefits of the accuracy of the growth models in describing the epidemic, along with a reasonable epidemiological interpretation of their free parameters. An attempt in this direction was presented by Wang [16] , where an approximate map between the Richards growth model [16] and the accumulated number of cases of a SIR model was proposed. The two free parameters of the Richards model were expressed as a function of the epidemiological based parameters of the SIR model. Here we improve on this analysis in two ways: (i) we extend the SIR model to a SIRD model by incorporating the deceased compartment, which is then used as the basis for the map onto the Richards model; (ii) the parameters of the SIRD model are allowed to have a time dependence, which is crucial to gain some efficacy in describing realistic cumulative epidemic curves of COVID-19. is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint

It is in general very hard to estimate the actual number of infected people within a given population, simply because a large proportion of infections go undetected. This happens largely because many carriers of the coronavirus are either asymptomatic or develop only mild symptoms, which in turn makes the number of confirmed cases for COVID-19 a poor proxy for the actual number of infections. This issue is well known in the literature and referred to as the "under-reporting problem" [17, 18] . With this in mind, and in the absence of more reliable estimates for the number of infected cases, we shall here focus our analysis on the fatality curves, defined as the cumulative number of deaths as a function of time.

In the present study we considered the mortality data of COVID-19 from the following countries: Italy, Germany, Sweden, Netherlands, Cuba, and Japan. The data used here were obtained from the database made publicly available by the Johns Hopkins University [19] , which lists in automated fashion the number of the confirmed cases and deaths attributed to COVID-19 per country. We have used data up to July 30, 2020.

The time evolution of the number of cases/deaths in an epidemy can be modelled by means of the Richards model (RM), defined by the following ordinary differential equation [20, 21, 22] :

where C(t) is the cumulative number of cases/deaths at time t, r is the growth rate at the early stage, K is the final epidemic size, and the parameter α measures the asymmetry with respect to the s-shaped curve of the standard logistic model, which is recovered for α = 1. In the present paper we shall apply the RM to the fatality curves of COVID-19, so that C(t) will always represent the cumulative numbers of deaths at time t, where t will be counted in days from the first death. Equation (3.1) must be supplemented with a boundary condition, which can be either the initial time, t = 0, or the inflection point, t = t c , defined by the condition

A direct integration of (3.1) yields the following explicit formula:

which will be the basis of our analysis.

. CC-BY-NC-ND 4.0 International license It is made available under a perpetuity.

is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint

We start by recalling the standard Susceptible (S)-Infected (I)-Recovered (R)-Deceased (D) epidemiological model

where S(t), I(t), R(t), and D(t) are the number of individuals at time t in the classes of susceptible, infected, recovered, and deceased respectively; whereas N is the total number of individuals in the population. i.e., N = S(t)+I(t)+R(t)+D(t).

The initial values are chosen to be S(0) = s 0 , I(0) = i 0 , with s 0 + i 0 = N , and R(0) = 0 = D(0). The parameters γ 1 and γ 2 are the rates at which infected individual becomes recovered or deceased, respectively. We then consider the following modified SIRD model, where in (3.3) and (3.4) we replace N with only the partial population in the S and I compartments, which takes into account the fact that the recovered (assuming they become immune) and the deceased cannot contribute to the transmission. We thus find

(3.10)

A fundamental quantity in epidemiology is the basic reproductive ratio, R 0 , which is defined as the expected number of secondary infections caused by an infected individual during the period she (or he) is infectious in a population consisting solely of susceptible individuals. In this model, R 0 can be calculated using the next generation method [23, 24] and is given by

Next, we define y(t) = S(t) + I(t) and divide (3.8) by (3.7) to obtain

. CC-BY-NC-ND 4.0 International license It is made available under a perpetuity.

is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint

Integrating both sides of (3.12), and inserting the result into (3.7), yields a growth equation of the Richards type:

where α = 1 − 1/R 0 and L = (i 0 + s 0 ) 1/α s 1−1/α 0 . We now seek to approximate the curve of accumulated death D(t), obtained from the SIRD model, with the Richards function C(t), as defined in (3.2) . To this end, we first impose the boundary conditions K = D(∞) and t c = t i , where t i is the inflection point of D(t). By definitionD(t i ) = 0, which implies from (3.10) thaṫ I(t i ) = 0 and thus

(3.14)

Furthermore, we require that at t = t i both C(t) and its derivativeĊ(t) respectively match D(t) andḊ(t), thus

Using the conditionİ(t i ) = 0 in the SIRD equations, we find

Using equations (3.15) and (3.16), we finally obtain the connection between the parameters (r, α) of the RM and the parameters (β, γ 1 , γ 2 ) of the SIRD model

, (3.20) which are the central equations of this paper. We can estimate the precision of the above 'map' between the RM and the SIRD model via the relative error function:

We have verified numerically that

. CC-BY-NC-ND 4.0 International license It is made available under a perpetuity.

is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint where ε is typically of order 0.1. A typical example of the agreement between the SIRD model, for a given set of parameters (β, γ 1 , γ 2 ), and the RM with the parameters obtained from the map described by (3.19) and (3.20) , is illustrated in Fig. 1. In Fig. 2 we show the simple monotonic dependence of the Richards parameters (r, α) on the parameter β of the SIRD model, for the biologically relevant interval 0 ≤ r, α ≤ 1. We also show, for comparison, the behavior of the basic reproduction number R 0 .

The SIRD model with constant parameters proved to be insufficient to accommodate properly the human intervention biased dynamics of the COVID-19 epidemics. The simplest solution to this problem is to allow the epidemiological parameter β to change in time according to the simple exponential decay function [4] 

where τ 0 is the starting time of the intervention and τ 1 is the average duration of interventions. Here β 0 is the initial transmission rate of the pathogen and the product β 0 β 1 represents the transmission rate at the end of the epidemic. Remarkably, the central map equations, (3.19) and (3.20) , are still valid, although t i is no longer given by (3.14) and should be determined from the maximum of the curve I(t) obtained from the numerical solution of the SIRD equations, with the parameter β replaced by the function β(t).

. CC-BY-NC-ND 4.0 International license It is made available under a perpetuity.

is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint Figura 2: Behavior of the Richards parameters (r, α) and R 0 as a function of the parameter β of the SIRD model.

In Fig. 3 we demonstrate some applications of the SIRD-RM map by showing the cumulative number of deaths (red circles) attributed to COVID-19 for the following countries: Italy, Germany, Netherlands, Sweden, Japan and Cuba. In all figures shown, the continuous (black) curve is the numerical fit to the empirical data, as produced by the SIRD model with the time-dependent parameter β(t) given in (3.23) , and the dashed (bright green) curve is the corresponding theoretical curve predicted by the RM, with the parameters as obtained from the map (3.19) and (3.20) . The statistical fits were performed using the Levenberg-Marquardt algorithm [25], as implemented by the lmfit Python package [26], to solve the corresponding non-linear least square optimization problem. In other words, the lmfit package was applied to each empirical dataset to determine the parameters (β 0 , β 1 , γ 1 , γ 2 , τ 0 , τ 1 ) of the SIRD model described in Secs. 3.2 and 3.3.

One can see from Fig. 3 that the agreement between the RM and the SIRD model is very good in all cases considered, which satisfactorily validates the map between these two models. This result thus shows, quite convincingly, that the parameters of the Richards model do bear a direct relationship to epidemiological parameters, as represented, say, in compartmental models of the SIRD type. Although the interpretation of the Richards parameters (r, α) are less obvious, in that they involve a nonlinear relation with the probability rates used in compartmental models, these parameters should nonetheless be regarded as bonafide epidemiological parameters. Furthermore, it is important to emphasize the flexibility of the RM: this model, which has only two time-independent parameters, is equivalent (in the sense of the is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint words, the two constant parameters of the RM are sufficient to characterize, to a rather good extent, the entire evolution of the COVID-19 epidemic in a given location.

It is worth pointing out that the discovery of power-law behaviors in the earlygrowth regime as well as in the saturation phase of the accumulated death curves, both of which are well described by the beta logistic model (BLM) [15] , brings about the challenge to accommodate power laws into a compartmental model. A preliminary analysis [15] shows that substantial modifications in the SIRD equations may be required to achieve power law-behavior in the short-time and long-time regimes of the epidemic curves. The possibility of a map between the BLM and a modified SIRD model with time-dependent parameters is currently under investigation.

The present paper provides a map between a SIRD model with time dependent parameters and the Richards growth model. We illustrated the use of this map by fitting the fatality curves of the COVID-19 epidemics data for Italy, Germany, Sweden, Netherlands, Cuba and Japan. The results presented here are relevant in that they showcase the fact that phenomenological growth models, such as the Richards model, are valid epidemiological models not only because they can successfully describe the empirical data but also because they capture, in an effective way, the underlying dynamics of an infectious disease. In this sense, the free parameters of growth models acquire a biological meaning to the extent that they can be put in correspondence (albeit not a simple one) with parameters of compartmental model, which have a more direct epidemiological interpretation.

[25] J. J. Moré, "The Levenberg-Marquardt algorithm: implementation and theory,"

in Numerical analysis, pp. 105-116, Springer, 1978.

[26] M. Newville, T. Stensitzki, D. Allen, and A. Ingargiola, "Non-linear leastsquares minimization and curve-fitting for Python," Chicago, IL, 2015.

. CC-BY-NC-ND 4.0 International license It is made available under a perpetuity.

is the author/funder, who has granted medRxiv a license to display the preprint in (which was not certified by peer review) preprint

The copyright holder for this this version posted August 6, 2020. . https://doi.org/10.1101/2020.08.04.20168120 doi: medRxiv preprint

Director-General's opening remarks at the media briefing on COVID-19 -30

Director-General's opening remarks at the media briefing on COVID-19 -30

Containing papers of a mathematical and physical character

Chinese and Italian COVID-19 outbreaks can be correctly described by a modified SIRD model

Power laws in superspreading events: Evidence from coronavirus outbreaks and implications for SIR models

Epidemiological model with anomalous kineticsthe Covid-19 pandemics

What will be the economic impact of covid-19 in the us? rough estimates of disease scenarios

Estimation of COVID-19 dynamics "on a back-of-envelope": Does the simplest SIR model provide quantitative parameters and predictions?

Modified SEIR and AI prediction of the epidemics trend of COVID-19 in China under public health interventions

Modeling the epidemic dynamics and control of COVID-19 outbreak in China

A SIR model assumption for the spread of COVID-19 in different communities

Generalized logistic growth modeling of the COVID-19 outbreak in 29 provinces in China and in the rest of the world

Modelling fatality curves of COVID-19 and the effectiveness of intervention strategies

Dynamics and future of SARS-CoV-2 in the human host

Complexity signatures in the COVID-19 epidemic: power law behaviour in the saturation regime of fatality curves

Richards model revisited: Validation by and application to infection dynamics

Age-structured estimation of COVID-19 ICU demand from low quality data

Analysis of COVID-19 under-reporting in Brazil

Coronavirus COVID-19 Global Cases by the Center for Systems Science and Engineering

A flexible growth function for empirical use

Richards model revisited: Validation by and application to infection dynamics

Richards model: a simple procedure for real-time prediction of outbreak severity

On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations

The construction of nextgeneration matrices for compartmental epidemic models

This work was partially supported by the National Council for Scientific and Technological Development (