Renormalization Group Approach to Pandemics as a Time-Dependent SIR Model

1 Introduction

Epidemic dynamics is often described in terms of a simple model introduced long time ago in [1]. Here, the affected population is described in terms of compartmentalised sub-populations that have different roles in the dynamics. Then, differential equations are designed to describe the time evolution of the various compartments. The sub-populations can be chosen to represent (S)usceptible, (I)nfected and (R)ecovered individuals (SIR model), obeying the following differential equations:

$\begin{array}{l}\frac{dS}{dt}=-\gamma S\frac{I}{N},\\ \frac{dI}{dt}=\gamma S\frac{I}{N}-\epsilon I,\\ \frac{dR}{dt}=\epsilon I,\end{array}(1.1)$

with the conservation law

$S\left(t\right)+I\left(t\right)+R\left(t\right)=N.(1.2)$

The system depends on three parameters, namely γ, ε and N. Due to the conservation law (2), only two equations are independent, so that one can drop the equation for S.

The cumulative number of infected, $\tilde{I}\left(t\right)$, that we are interested in, is related to the above sub-populations as

$\tilde{I}\left(t\right)=I\left(t\right)+R\left(t\right).(1.3)$

We can therefore re-write the two independent SIR equations as

$\frac{d\tilde{I}\left(t\right)}{dt}=\gamma \left(\tilde{I}\left(t\right)-R\left(t\right)\right)\left(1-\frac{\tilde{I}\left(t\right)}{N}\right),(1.4)$

$\frac{dR\left(t\right)}{dt}=\epsilon \left(\tilde{I}\left(t\right)-R\left(t\right)\right).(1.5)$

Empirical modifications of the basic SIR model exist and range from including new sub-populations to generalise the coefficients γ, ε to be time-dependent in order to better reproduce the observed data (see Refs. [2–5] for recent examples).

Recently the epidemic Renormalization Group approach (eRG) to pandemics, inspired by particle physics methodologies, was put forward in [6]. The approach builds on scale transformations and invariances as introduced by Wilson in [7, 8], applied in this case to the evolution of social systems. The method has been further explored in [9], where it was demonstrated the eRG effectiveness when describing how the pandemic spreads across different regions of the world. This allowed to successfully simulate the onset of the second wave pandemic in Europe [10]. The eRG framework provided useful also to determine the impact and level of social distancing when combined with mining the Google and Apple mobility data [11]. We refer the reader to the review in [12] for a clear and comparative discussion of several approaches used in order to model the COVID-19 pandemic.

The goal of the present work is to further extend the eRG framework to properly take into account the number of recovered cases so that a better understanding of the reproduction number can also be achieved. We will start, first, by providing a map between the original eRG model and certain modified SIR models, establishing in this way and for the first time an exact connection between two powerful approaches to pandemics. We will finally test the framework via COVID-19 data.

1.1 Reviewing the Epidemic Renormalization Group

In the epidemic Renormalization group (eRG) approach [6], rather than the number of cases, it is convenient to discuss its logarithm, which is a more slowly varying function

$\alpha \left(t\right)=\text{ln}\tilde{I}\left(t\right),(1.6)$

where $\text{ln}$ indicates the natural logarithm. The derivative of α with respect to time provides a new quantity that we interpret as the beta-function of an underlying microscopic model. In statistical and high energy physics, the latter governs the time (inverse energy) dependence of the interaction strength among fundamental particles. Here it regulates infectious interactions.

More specifically, as the Renormalization group equations in high energy physics are expressed in terms of derivatives with respect to the energy μ, it is natural to identify the time as $t/t_0=-\text{ln}\left(\mu /{\mu }_0\right)$, where $t_0$ and ${\mu }_0$ are respectively a reference time and energy scale. The latter is introduced to ease the reading and notation for the physics community. We choose $t_0$ to be one week so that time is measured in weeks, and will drop it in the following. Thus, the dictionary between the eRG equation for the epidemic strength α and the high-energy physics analog is

$\beta \left(\alpha \right)=\frac{d\alpha }{d\ln \left(\mu /{\mu }_0\right)}=-\frac{d\alpha }{dt}.(1.7)$

It has been shown in [6] that α captures the essential information about the infected population within a sufficiently isolated region of the world. The pandemic beta function can be parametrised as

$-\beta \left(\alpha \right)=\frac{d\alpha }{dt}=\tilde{\gamma }\alpha {\left(1-\frac{\alpha }{a}\right)}^n,(1.8)$

whose solution, for $n=1$, is a familiar logistic-like function

$\alpha \left(t\right)=\frac{a{e}^{\tilde{\gamma }t}}{b+e^{\tilde{\gamma }t}}.(1.9)$

The dynamics encoded in Eq. (1.8) is that of a system that flows from an UV fixed point at $t=-\infty$ where $\alpha =0$ to an IR fixed point where $\alpha =a$. The latter value encodes the cumulative number of infected cases $P=\text{exp}\left(a\right)$ in the region under study at the end of the epidemic wave. The coefficient $\tilde{\gamma }$ is the diffusion slope, while b shifts the entire epidemic curve by a given amount of time. Further details, including what parameter influences the flattening of the curve and location of the inflection point and its properties can be found in [6].

1.2 Connecting Epidemic Renormalization Group With SIR While Extending It

To start connecting with compartmental models we rewrite Eq. (1.8) as

$\frac{d\tilde{I}}{dt}=\tilde{\gamma }\tilde{I}\ln \tilde{I}\left(1-\frac{\ln \tilde{I}}{\ln P}\right),(1.10)$

whose solution, with the initial condition $\ln {\tilde{I}}_0=\ln \tilde{I}\left(0\right)=\frac{\ln P}{b+1}$, is a logistic function written as

$\tilde{I}\left(t\right)=\text{exp}\left(a\frac{1}{b{e}^{-\tilde{\gamma }t}+1}\right).(1.11)$

In the original eRG framework the number of recovered individuals were not explicitly taken into account. This is, however, straightforward to implement by introducing an equation for $dR/dt$ and imposing a conservation law equivalent to the one for the SIR model. A minimal choice compatible with the conservation law is

$\begin{array}{l}\frac{dR}{dt}=\epsilon I,\\ \frac{dI}{dt}=\tilde{\gamma }\left(I+R\right)\text{ln}\left(I+R\right)\left(1-\frac{\text{ln}\left(I+R\right)}{\ln P}\right)-\epsilon I,\\ \frac{dS}{dt}=-\tilde{\gamma }\left(I+R\right)\text{ln}\left(I+R\right)\left(1-\frac{\text{ln}\left(I+R\right)}{\ln P}\right),\end{array}(1.12)$

where the parameters are $\tilde{\gamma }$, P and ε. At fixed N, $\tilde{\gamma }$, P and for any value of ε, the SIR model in (1.1) and the eRG systems of equations match if we allow γ to be the following time-dependent function

$\gamma \left(t\right)=\tilde{\gamma }\left(I\left(t\right)+R\left(t\right)\right)\text{ln}\left(I\left(t\right)+R\left(t\right)\right)\left(1-\frac{\text{ln}\left(I\left(t\right)+R\left(t\right)\right)}{\ln P}\right)\frac{N}{I\left(t\right)S\left(t\right)}.(1.13)$

As we shall see this is a welcome feature. To better appreciate the mapping we show in Figure 1 the time-dependent γ parameter for a hypothetical case with $N=7$ millions, $P=300$ thousands, $\epsilon =1$ and $\tilde{\gamma }=0.7$ with initial conditions $I\left(0\right)=3$, $R\left(0\right)=0$ and $S\left(0\right)=N-I\left(0\right)$.

The result is a smooth function that peaks at short times and then plateaus to a fraction of $\tilde{\gamma }$. In other words the eRG naturally encodes a rapid diffusion of the disease in the initial states of the epidemic and the slow down at later times.

In terms of $\tilde{I}$ one can more compactly write the system of equations in (1.12) as:

$\begin{array}{l}\frac{dR}{dt}=\epsilon \left(\tilde{I}-R\right),\\ \frac{d\tilde{I}}{dt}=\tilde{\gamma }\tilde{I}\ln \tilde{I}\left(1-\frac{\ln \tilde{I}}{\ln P}\right),\\ \frac{dS}{dt}=-\frac{d\tilde{I}}{dt}.\end{array}(1.14)$

FIGURE 1

FIGURE 1. Typical form of time-dependent γ matching the SIR system to the eRG one.

1.3 Reproduction Number

An important quantity for pandemics is the reproduction number related to the expected average number of infected cases due to one case. In the time-dependent generalised SIR model it is identified as:

$R_0\left(t\right)=\frac{\gamma \left(t\right)}{\epsilon \left(t\right)}.(1.15)$

To extract $R_0$ from data it is useful to recast it as:

$R_0=\frac{\gamma }{\epsilon }=\frac{\frac{dI}{dt}+\frac{dR}{dt}}{\frac{dR}{dt}}\frac{N}{S}=\frac{\frac{d\tilde{I}}{dt}}{\frac{dR}{dt}}\frac{N}{S}.(1.16)$

The result holds for γ and ε generic functions of time. Here we also generalise the time-dependence of ε to be:

$\epsilon \left(t\right)=A\left[1-c\cdot e^{-\frac{1}{2}{\left(\frac{t-t_0}{W}\right)}^2}\right].(1.17)$

The choice above has been empirically devised to best describe the data within the current approach. As it is clear from its form this function has a dip at $t_0$ (possibly correlated with the peak of the newly infected cases) of width W and depth $c\cdot A$ with A the asymptotic value for $t\to \pm \infty$.

As we shall see the shape allows for a substantial increase of $R_0$ near the peak of the newly infected cases that could be due to a number of factors including possible health-system stress around this period.

2 Testing the Framework

As a timely application we consider the COVID-19 pandemic. Here the factor $N/S$ in (1.16) can be neglected as the number of total infected is at most of $\text{O}\left(1\text{\%}\right)$ of the total susceptible population and the ratio is therefore very close to unity. The reproduction number can hence be estimated as:

$R_0\left(t\right)=\frac{\text{newly}\text{infected}\left(\text{t}\right)}{\text{newly}\text{recovered}\left(\text{t}\right)}.(2.1)$

The values for the numerator and denominator for different regions of the world can be obtained from several sources such as the World Health Organization (WHO) and Worldometers.

As testbed scenarios we consider four benchmark cases, namely Denmark, Germany, France and Italy. These countries adopted different degrees of containment measures. We find convenient to bin the data in weeks to smooth out daily fluctuations. We associate an error to both newly infected and newly recovered given by the square root of their values.

Procedurally we first fit the function $\alpha \left(t\right)$ to determine a, b and $\tilde{\gamma }$ following [6, 9]. We then solve with these, as input, the system of equations (1.12) for $\text{O}\left({10}^4\right)$ different choices, within reasonable ranges, of the parameters A, c, $t_0$ and W entering the definition of $\epsilon \left(t\right)$, Eq. (1.17). The optimal choice of such parameters is finally obtained by performing a ${\chi }^2$ minimization to the data related to the recovered cases. Combining the results with $\gamma \left(t\right)$ we compare $R_0\left(t\right)$ with the actual data.

For each country we show the data and the model results by grouping together five graphs in a single figure. The different panels represent $\tilde{I}\left(t\right)$, $R\left(t\right)$, $R_0\left(t\right)$, $\gamma \left(t\right)$ and $\epsilon \left(t\right)$, all as function of the week number. Additionally the data will be reported starting some time after the outbreak. The reason being that the values of the number of recovered cases at early times is too small to be reliable and begins to be sizeable only few weeks after the outbreak. For our predictions $\epsilon \left(t\right)$, $\gamma \left(t\right)$ and $R_0\left(t\right)$ we show bands limiting the 90% confidence level. Those are obtained shifting the data for the number of recovered cases by 1.65 standard deviations. The fitting errors for $\tilde{\gamma }$, a and b from the method in [6] can be neglected given that these parameters are highly constrained by the data.

In general we find good agreement between the data and the model with the exception of the increase in the number of infected cases occurring in the last weeks for some of the countries. Those are non-smooth events resulting from an abrupt change in the social distancing measures or from new and previously un-accounted disease hotspots. As such those events cannot be predicted by smooth models.

2.1 Denmark

The data and the model results for Denmark are shown in Figure 2. For the time-dependence of $R_0$, given in panel 2c, we observe that the model captures the variation of the data for over 10 weeks. Overall we find that the eRG model provides a reasonable description of the time dependence of the reproduction number. We observe that the recovery rate $\epsilon \left(t\right)$ grows with time by a factor of five. There could be several factors contributing to this growth, one being a better trained health system.

FIGURE 2

FIGURE 2. Time dependence of $\tilde{I}\left(t\right)$, $R\left(t\right)$, $R_0\left(t\right)$, $\gamma \left(t\right)$ and $\epsilon \left(t\right)$ within the eRG SIR model optimised to describe the Danish data for COVID-19. Solid lines are the model (shown as a 95% confidence level band for the predictions) and the dots with error-bars are the data. For the chosen confidence intervals the errors for the fit functions $\tilde{I}$ and R are not visible in this and the following figures.

2.2 Germany

For Germany the analysis is summarised in Figure 3. The overall trends are similar to the Danish case including the temporal trend of the recovery rate $\epsilon \left(t\right)$.

FIGURE 3

FIGURE 3. Time dependence of $\tilde{I}\left(t\right)$, $R\left(t\right)$, $R_0\left(t\right)$, $\gamma \left(t\right)$ and $\epsilon \left(t\right)$ within the eRG SIR model optimised to describe the German data for COVID-19. Solid lines are the model (shown as a 95% confidence level band for the predictions) and the dots with error-bars are the data.

2.3 Italy

The analysis for Italy is shown in Figure 4. We observe rather large values of $R_0$ compared to Denmark at early times and a factor of two with respect to Germany. We also observe that a good fit is obtained for $\gamma \left(t\right)$ approaching very small values at large times. This is different from Germany and Denmark, suggesting strong distancing measures being adopted by the Italian government. This seems to be further followed by a smaller value of the recovery rate, roughly about a fourth. However this last comparison is biased by the fact that the number of deaths in Italy is about 15% of the number of infected cases while in Germany and Denmark it is below 5% suggesting that a more accurate description at large times would require introducing also a compartment accounting for the deaths.

FIGURE 4

FIGURE 4. Time dependence of $\tilde{I}\left(t\right)$, $R\left(t\right)$, $R_0\left(t\right)$, $\gamma \left(t\right)$ and $\epsilon \left(t\right)$ within the eRG SIR model optimised to describe the Italian data for COVID-19. Solid lines are the model (shown as a 95% confidence level band for the predictions) and the dots with error-bars are the data.

2.4 France

The results for the French case can be found in Figure 5 with the overall picture similar to the Italian one. The striking difference compared to the other countries is that the recovery rate decreases at late times. Given that the deaths in France amount to about 20% of the total infected, such difference indicates that a more complete model (including the deaths compartment) is needed.

FIGURE 5

FIGURE 5. Time dependence of $\tilde{I}\left(t\right)$, $R\left(t\right)$, $R_0\left(t\right)$, $\gamma \left(t\right)$ and $\epsilon \left(t\right)$ within the eRG SIR model optimised to describe the French data for COVID-19. Solid lines are the model (shown as a 95% confidence level band for the predictions) and the dots with error-bars are the data.

3 Conclusions and Outlook

We generalised the epidemic Renormalization Group framework to take into account the recovered cases and to be able to determine the time dependence of the reproduction number. At the same time we show that the eRG framework can be embedded into a SIR model with time-dependent coefficients. Interestingly the resulting infection rate $\gamma \left(t\right)$ is a smooth curve with a maximum at early times while rapidly plateauing at large times. This is a welcome behaviour since it encodes the slow down in the spreading of the disease at large times coming, for example, from social distancing.

We then move to confront the model to actual data by considering the spread of COVID-19 in the following countries: Denmark, Germany, Italy and France. We show that the overall approach works rather well in reproducing the data. Nevertheless the interpretation for the recovery rate is natural for Denmark and Germany while it requires to add the deaths compartment for Italy and France. The reason being that for the first two countries the number of deaths is below 5% of the number of infected cases while it is above 15% for France and Italy¹. We therefore expect that this compartment will be relevant to include for countries with similar number of deaths. The extension of the eRG to include also the death compartment is a natural next step. Upon introducing thD-compartment, the system of differential equations in 1.1 (and similarly for Eq. 1.12) extends to

$\begin{array}{l}\frac{dS}{dt}=-\gamma S\frac{I}{N},\\ \frac{dI}{dt}=\gamma S\frac{I}{N}-\left(1-\alpha \right)\epsilon I-\alpha \rho I,\\ \frac{dR}{dt}=\left(1-\alpha \right)\epsilon I,\\ \frac{dD}{dt}=\alpha \rho I.\end{array}(3.1)$

with $I\left(t\right)+R\left(t\right)+S\left(t\right)+D\left(t\right)$. becoming the new conserved quantity. Notice that two new parameters had to be introduced, i.e. the fatality ratio α and the death rate ρ. In the case of COVID-19 the data suggest that these parameters have an important time dependence, as it is manifest from Figure 6.

FIGURE 6

FIGURE 6. The average fatality ratio as a function of the cumulative number of infected per million inhabitants for the countries discussed here.

Here the running average fatality ratio (cumulative number of deaths over cumulative number of infected) is displayed, as a function of the number of infected per million inhabitants, for the countries considered here.

We are now investigating different parameterizations of these curves, however, since our main focus here is to establish an exact connection between two well known frameworks describing social dynamics (the renormalization group and compartmental models) we have restricted our attention to the simplest compartmental model, i.e., the SIR model. To our knowledge such an exact matching could not be found in previous literature. Another interesting approach could be to include the death compartment in the recovered one as put forward in [13].

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: https://www.who.int/csr/sars/country/en/ World Health Organization (WHO) and https://www.worldometers.info/coronavirus/Worldometers.

Author Contributions

All authors listed have made a substantial, direct, and intellectual contribution to the work and approved it for publication.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The reviewer KT declared a past co-authorship with one of the authors FS to the handling editor.

References

1. Kermack WO, McKendrick AG, A contribution to the mathematical theory of epidemics, Proc R Soc A (1927) 115(772):700–21. doi:10.1007/BF02464423

CrossRef Full Text | Google Scholar

2. Zheng M, Zhao M, Min B, Liu Z. Synchronized and mixed outbreaks of coupled recurrent epidemics. Sci Rep (2017) 7:2424. doi:10.1038/s41598-017-02661-9

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Zheng M, Wang C, Zhou J, Guan S, Zou Y, et al. Non-periodic outbreaks of recurrent epidemics and its network modelling. Sci Rep (2015) 5:16010. doi:10.1038/srep16010

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Li R, Pei S, Chen B, Song Y, Zhang T, Yang W, et al. Substantial undocumented infection facilitates the rapid dissemination of novel coronavirus (SARS-CoV-2). Science (2020) 368(6490):489–93. doi:10.1126/science.abb3221

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Zhang X, Ma R, Wang L, Predicting turning point, duration and attack rate of COVID-19 outbreaks in major Western countries. Chaos Solit Fractals (2020) 140:110130. doi:10.1016/j.chaos.2020.109829

CrossRef Full Text | Google Scholar

6. Della Morte M, Orlando D, Sannino F, Renormalization group approach to pandemics: the COVID-19 case. Front Physiol (2020) 8:144. doi:10.3389/fphy.2020.00144

CrossRef Full Text | Google Scholar

7. Wilson K, Renormalization group and critical phenomena. I. Renormalization group and the Kadanoff scaling picture. Phys Rev B (1971) 4:3174–83. doi:10.1103/PhysRevB.4.3174

CrossRef Full Text | Google Scholar

8. Wilson K, Renormalization group and critical phenomena. II. Phase-space cell analysis of critical behavior. Phys Rev B (1971) 4:3184–205. doi:10.1103/PhysRevB.4.3184

CrossRef Full Text | Google Scholar

9. Cacciapaglia G, Sannino F, (2005) [arXiv:2005.04956 [physics.soc-ph]]. Accepted for publication in Nature Scientific Reports.

Google Scholar

10. Cacciapaglia G, Cot C, Sannino F, (2007) [arXiv:2007.13100] [physics.soc-ph]].

Google Scholar

11. Cacciapaglia G, Cot C, Sannino F, (2008) [arXiv:2008.02117] [physics.soc-ph]].

Google Scholar

12. Smit AJ, Fitchett JM, Engelbrecht FA, Scholes RJ, Dzhivhuho G, Sweijd NA, Winter is coming: a southern hemisphere perspective of the environmental drivers of SARS-CoV-2 and the potential seasonality of COVID-19, Int J Environ Res Publ Health (2020) 17:5634. doi:10.3390/ijerph17165634

CrossRef Full Text | Google Scholar

13. Scala A, Flori A, Spelta A, Brugnoli E, Cinelli M, Quattrociocchi W, et al. Time, space and social interactions: exit mechanisms for the Covid-19 epidemics Sci Rep (2020) 10:13764. doi:10.1038/s41598-020-70631-9

PubMed Abstract | CrossRef Full Text | Google Scholar