key: cord-0800105-sad09fby
authors: Duclos, T. G.; Reichert, T. A.
title: The Missing Link: A Closed Form Solution to the Kermack and McKendrick Model Equations
date: 2021-03-05
journal: nan
DOI: 10.1101/2021.03.02.21252781
sha: d0e48c0e19eee92ffdc2a612d6c4ebc97f9ae471
doc_id: 800105
cord_uid: sad09fby

Susceptible infectious recovered (SIR) models are widely used for estimating the dynamics of epidemics. Such models project that containment measures flatten the curve, i.e., reduce but delay the peak in daily infections, cause a longer epidemic, and increase the death toll. These projections have entered common understanding; individuals and governments often advocate lifting containment measures such as social distancing to shift the peak forward, limit societal and economic disruption, and reduce mortality. It was, then, an extraordinary surprise to discover that COVID-19 pandemic data exhibit phenomenology opposite to the projections of SIR models. With the knowledge that the commonly used SIR equations only approximate the original equations developed to describe epidemics, we identified a closed form solution to the original epidemic equations. Unlike the commonly used approximations, the closed form solution replicates the observed phenomenology and quantitates pandemic dynamics with simple analytical tools for policy makers. The complete solution is validated using independently measured mobility data and accurately predicts COVID19 case numbers in multiple countries.

The susceptible-infectious-recovered (SIR) models (1) , widely used in predictions of epidemic spread, are said to predict that social distancing will "flatten the curve", i.e., reduce and delay the peak of new daily cases. These models project that the peak in daily cases will be reduced by a modest fraction, but the interval within which case numbers remain large lengthens in approximately inverse proportion. With this image highlighted in the popular media, concerns about the economic devastation projected to be associated with such measures applied in response to COVID-19 have caused many individuals and governments to advocate that social distancing measures to be lifted early (2) to shift the peak forward, to end the economic, educational, and social disruption sooner, and to reduce mortality.

As the COVID-19 pandemic swept across the world in the winter and spring of 2020, different countries applied diverse mitigation measures (3) (4) (5) in their attempts to control the spread of the virus. Although not deliberate, these countries have conducted "natural experiments" on the effectiveness of differing levels of containment measures. Their case data provides an historically unique opportunity to compare different models' projections of the effect of diverse containment measures against real-world data.

In this manuscript, we first compare the trends in COVID-19 case data from countries that applied differing levels of containment to the trends predicted by an SIR model. We show that SIR models do not replicate even the basic trends in the country data. This discovery prompted us to develop a closed form solution to the original equations developed by Kermack and McKendrick (1) to which the SIR models are an approximation. We verify the solution by accurately predicting time series of case data for a sample of countries and by correlating the . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint model to independent data. Finally, we show how to use the solution to manage epidemics with quantitative tools.

In 1927, Kermack and McKendrick (1) developed the following complex set of integrodifferential equations for modelling epidemics:

(1, S1-5)

(2, S1-6) as "the rate of removal" of the infected population, and ܰ ௗ is defined as the initial population . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. Although Kermack and McKendrick attempted to find an analytical solution for the CSIR equations, they, and many other researchers since (6), were only successful in finding solutions to approximations of their equations. One such set of approximate equations proposed by the original authors is shown here:

(6, S1-1)

, and (8, S1-3) ܰ ൌ ܵ ‫ܫ‬ ܴ , (9, S1-4)

where β = rate of contact and transmission; γ = rate of recoveries = 1/t r ; and t r = time of infectiousness.

These equations (equations 6-9) are the well-known "SIR" model equations, which we will call here the Approximate SIR (ASIR) model equations. These equations can be derived from . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. , is the criterion for when the epidemic is beginning to end, and therefore, a useful model should be well-behaved at and in the vicinity of this condition.

Despite the mathematical problems with equations 6-9, the true test of a model, even an approximate one, is the extent to which it predicts real data trends. It is fortuitous, then, that the progression of the COVID-19 pandemic has been well documented in multiple countries which took different paths while attempting to contain the spread of the virus. This dataset affords us the opportunity to test the veracity of the trends predicted by ASIR model against actual data.

As a first step, we used equations 6-9 to project trends in the COVID-19 pandemic. In Figure 1A , we see that the ASIR model projects that the time at which the total number of cases levels off is delayed as social distancing increases (represented by decreasing β ), while the ultimate total case number remains similar regardless of the value of β . In this simulation, increasing social distancing also delays and broadens the peak of new infections ( Figure 1B ).

In a second step, we compared the ASIR projected trends to case data (8) for the COVID-19

pandemic in Sweden and New Zealand ( Figure 1C , D), and for South Korea and Italy ( Figure   1E , F). These pairs of countries have comparable population densities, but implemented . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint mitigation measures with different timings and intensities (3) (4) (5) 9) . In particular, New Zealand and South Korea introduced stronger social distancing measures much earlier than Italy and Sweden.

In contrast to the ASIR model predictions, the country data in Figures 1C-F show that stronger social distancing measures are associated with an earlier and lower peak in new infections and an earlier levelling off at a lower total number of cases. Both trends demonstrate that the ASIR model ( Figure 1A , B) is not merely inaccurate but projects epidemic data to trend in the opposite direction to the reported data. Other authors (10), too, have noted that the peaks of cases in countries with stronger social distancing occur earlier than in countries with weaker social distancing.

The ASIR model fails the simplest test of model veracity: the projection of qualitative trends. McKendrick model equations. This solution was developed by first recasting the original equations into a new set of differential equations, and then solving these equations using basic . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint principles. The closed form solution incorporates measurable population interactions, and can be manipulated using straightforward mathematical techniques to reveal a number of relationships useful for managing an epidemic.

As a first step towards a solution, the CSIR equations (equations 1-5) were recast into the following set of differential equations:

(13, S1-29)

(14, S1-13) using these definitions for

(15, S1-24)

, and ‫ܣ‬ ஶ is the area that the population ܰ ஶ inhabits.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint As detailed in Supplement 1, equations 10-16 were developed by adopting a perspective from inside the epidemic rather than using an external view. Equations 10-16 are also shown to be mathematically identical to equations 1-5 in Supplement 1, but recast in a form for which a solution can be found.

In a second step, using basic principles (Supplement 1) and equations 10, 13 and 14, we developed the following differential equation to describe the new daily cases during an epidemic: S1-65) and, in Supplement 1, the following are shown to be the solutions to equations 10, 11 and 12:

and (20, S1-68)

(21, S1-69)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) 

(24, S1-83)

Equations 18-24 are a closed form solution to the original model proposed by Kermack and McKendrick and, as we will show in the next section, they accurately describe the evolution of an epidemic.

The qualitative behavior of the CSIR solution is plotted in Figure 1G , H. The CSIR solution predicts that with increased social distancing (higher K 1 ), the total number of daily cases is lower and attained earlier ( Figure 1G ). The peak in new cases is also lower and occurs earlier ( Figure   1H ). Unlike the ASIR model trends ( Figure 1A , B), the CSIR solution trends ( Figure 1G , H) are identical to those exhibited by the real-life COVID-19 data in Figure 1C -F. Whereas the ASIR model fails to predict the proper impact of social distancing, the CSIR solution predicts the trends seen in the real data and, in particular, that stronger social distancing produces an earlier peak of daily cases. We also note here that, relative to the native, unmanaged epidemic, for all containment measures (higher K 1 ), the peak is lower and occurs earlier. In other words, the curve never flattens; rather, it reaches a peak earlier and then falls in a steep decline.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) were then used to predict the balance of the epidemic data that followed.

Using this approach, the CSIR solution predicted the course of total cases ( Figure 3 ) with an R 2 > 0.97 for all six countries for the 45 days following the date containment measures were introduced (Table 1 ). In the case of Sweden, which did not introduce containment measures, March 23, 2020 was chosen as the starting point.

The CSIR solution also predicted daily new cases ( Figure 4 ) for these six countries for the 45 days following the date containment measures were introduced with an R 2 range of 0.29 to 0.90

As seen in the figure, the predicted peak of new cases was close to the observed peak for all countries.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) It is important to emphasize that the CSIR projections in Figures 3 and 4 are not fits to the fulllength of the data shown. Rather, a short portion of the epidemic data, with specified characteristics, was excised for the purposes of determining coefficients of the model; these were then used to project the rest of the curves before and after the time frame of the excised data.

In an additional demonstration of the CSIR solution's veracity, we tested the assumption that ‫ܭ‬ ் is a property of the disease; and, therefore, should be the same for each country. Equation S1-79 shows that the model parameters, expressed in a function, F(N(t)), should be linearly proportional to time with a constant of proportionality −K T . As illustrated in Figure 5 , the fit of equation S1-79 to the country data from Table 2 has an R 2 = 0.97 and a slope of −0.25 (which is equal to −K T ). This excellent correlation confirms that ‫ܭ‬ ் is the same for all countries; and, therefore, likely is a property of the disease.

A third illustration of the solution's veracity was created by correlating mobility data obtained over the modelled period for each of the countries. These data are available from Google (11) and are a measure of the difference between the amount of time people stayed at home (the "Residential" data set) during the period modelled and a baseline measured for 5 weeks starting January 3, 2020. Based on equation S1-46, the integral of this mobility data should correlate with the measured RCO. Figure 6 shows that for each country considered, the correlation was quite high.

These excellent results demonstrate that the CSIR solution, embodied by equations 18-24, correctly characterizes epidemic dynamics from multiple countries in a unified way, something that ASIR models simply cannot duplicate. They also demonstrate that the model developed by . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. If an epidemic is not controlled in its early stages, the tools explicated in Supplement 3 can be used to determine the state of the epidemic and to calculate quantitative actions that must be taken to control and end the epidemic. Strong intervention actions are still necessary and, using the RCO metric, the effectiveness of these actions can be quickly determined from the resultant . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint case trends. The epidemic state, tracked with the RCO metric, can be closely monitored and control measures adjusted according to the observed trends.

As an epidemic progresses, outbreaks and surges should be expected. These will occur if containment flags, more infectious variants emerge, or new cases are introduced from outside the region of focus. All will need to be controlled. As explained in Supplement 3, the start of an outbreak causes obvious changes in the behavior of the RCO metric. Diligently monitored, these RCO changes can be detected early enough for public health officials to react in a timely manner, to bring outbreaks under control, and to re-establish the planned course.

The data on the current COVID-19 pandemic are widely accessible from a variety of sources, updated daily. The unfolding panorama provides a test bed for models used to predict outcomes and the effects of various interventions. Because different countries have employed different containment strategies (3) (4) (5) 9) , the world is conducting an epidemiological experiment on a grand scale.

Current epidemiological models use the ASIR equations, which have been assumed to be reasonably accurate representations of the complete equations developed in 1927. They exist in many variants, both deterministic and stochastic; and their behavior is widely known. It is startling, then, that when the classic ASIR model is tested using the currently available data from the COVID-19 pandemic, it fails the most basic test for any model: it projects trends that are the opposite of those easily visible in the data from many countries.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint

Unfortunately, the ASIR model and its variants have been used to fashion guidelines for epidemic management. Tragically, the false notion that stronger social distancing lengthens the pandemic may well have caused country leaders, especially those most concerned with economic performance, to see social distancing as producing only a modest reduction in the horror of an epidemic peak while significantly prolonging economic disruption.

The CSIR solution, developed from basic mathematical principles, accurately projects the epidemic trends. It makes clear that short and sharp social distancing produces rapid truncation of epidemic upward trends, thereby shortening-not lengthening-the time needed to bring epidemics under control. Second, an indicator of the rate of change in epidemic dynamics (the RCO) allows direct observation of the effectiveness of intervention measures and provides policy makers with an opportunity to react before new outbreaks gain momentum.

Every country and economy can use the solution presented here to plan and implement the highest level of social distancing measures deemed sustainable to quickly reduce case numbers to levels at which case identification, contact tracing, testing, and isolation can be maintained, allowing a more rapid return to nearly normal social interactions while minimizing economic consequences. The ultimate insight from the model is one of hope: the path of an epidemic is not an uncontrollable force of nature; nor is epidemic control inevitably the road to economic ruin.

Rather, the afflicted population can, through their behavior, choose to control their destiny.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint dashed curves), equation 18 projects that total cases will rise to lower levels, and reach these levels earlier (G) Similarly, equation 17 projects that new daily cases will peak earlier and lower values with increasing containment (H).

The CSIR model trends in (G) and (H) are highly similar to those in the country data (C -F).

The ASIR model trends in (A) and (B) have completely different shapes; and vary with increasing containment in an opposite sense to that in the country data.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint any specific containment measures, so the model calibration was begun on April 1, when the slope of the rate of change operator (RCO)curve first became steady.). All dates are in 2020.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. March 23-May7. The CSIR model was calibrated using data from the date ranges listed in Table   re . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. Parameters from linear fit of rate of change operator (RCO) data in Figure 4 : K 1 , slope; K 2 , intercept; N ( t 0) , number of cases at time (t 0 ), first day of range used. All dates are in 2020.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint A solution to equations 1-4 is easily obtained using straightforward computational methods. Figure 1A and 1B depict one such solution computed using a simple Euler method to solve these differential equations with a time step of 1 day.

Because the assumption that β and ߛ are constant strongly limits the capability to describe epidemics, we call this form of the model the "approximate SIR" (ASIR) model. While not intended to model epidemic dynamics perfectly, equations 1-4 have always been assumed to reflect the general trends of an epidemic (1). As Figure 1 clearly demonstrates, however, these equations do not even qualitatively reproduce the dynamics of an epidemic.

Referring to the original equations (1), we can write the epidemic equations without any assumptions about the time-varying nature of any parameter, using our notation of S, I, and R from above as:

(S1-8, 4)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint to the current infections. The third step is to recognize that the values of the integrals in equations S1-5-7 are equal to zero and ‫ܤ‬ ሺ ߠ ሻ ൌ 1 whenever the epidemic starts. With these perspectives in mind, if equations S1-5-7 and S1-9 are first all multiplied by ‫ܣ‬ ௌ (the area that ܰ ௌ ሺ ‫ݐ‬ ሻ inhabits when ‫ݐ‬ ൌ 0 ), they can be rewritten as 

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. 

(S1-19)

Taking the limit as ∆ ‫ݐ‬ ՜ 0

, we obtain the following differential equations:

(S1-23) 

, will be dropped.

The perspective in the immediately preceding part of the analysis is that the epidemic can be considered to start over again at each instant in time. In this way, the susceptible population is . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint not fixed by initial conditions, but rather is the population that will eventually become infected in the future. Embedded in this concept is the assumption that, during each

, the susceptible population is always in contact with those people that have been previously infected or who will become infected, i.e., the epidemic remains contiguous. This is not restrictive when considering the initial stage of the epidemic; however, it plays a critical role in recognizing an outbreak that starts outside the population,

Equations S1-20-23 can also be developed from equations S1-5-7 using a different consideration. If we first define two functions

(S1-25, 16)

then we can rewrite equations S1-5-7 in the following form after multiplying by

and (S1-28, 12) (S1-29, 13) . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. equations S1-26-29 are also equivalent to equations S1-20-23. Both sets of equations are mathematically equivalent to the Kermack and McKendrick (1) system of equations, S1-5-9.

We in the solutions to equations S1-5-9 or S1-26-29.

The solution to the CSIR model

We now need to find a solution to equations S1-26-29. We begin by reiterating the definitions of = R(t) + I(t). The population, R(t), is also referred to as the recovered population because although they may still be ill, they are no longer infectious.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint

We assume that people in the recovered population remain immune after recovery and that the epidemic begins with the introduction of the infection by one individual at time t=0; therefore, I(0) = N(0) = 1. In addition, we assume that no new infections (other than the initial infection) are introduced from outside the region of interest during the epidemic.

The goal in solving equations S1-26-29 is to find an expression for N(t) and for the total number of people who will be infected during the epidemic, defined as ܰ ஶ , in terms of the parameters of disease transmission and the behavior of the population enduring the epidemic. We will first find an expression for N(t) and then find ܰ ஶ by taking the limit of N(t) as time increases. Along the way, we will also find expressions for I(t), R(t), and S(t).

The development of the solution begins by noting the following:

and (S1-32)

is a yet unknown function of time.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. 

can only change to exactly the extent that the quantity

changes. Thus, we can write the following expression:

is a function that modifies the fraction

in association with recoveries at time t. It is initially assumed to be a function of time.

Equation S1-34 is a simple differential equation in the variable

, whose solution is

.

(S1-35)

(S1-36)

Likewise, the solutions for R(t) are

, and (S1-37)

(S1-38)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. , (S1-39)

, and (S1-41)

(S1-42)

During an epidemic, if every contact made by a person within the population N(t) could result in an infection, then the disease progression can be represented by the following expression:

where ‫ܭ‬ is a constant representing the transmissibility of the disease.

consists of both infected and recovered individuals, we can rewrite equation S1-43

(S1-44)

Because the population R(t) cannot transmit the disease, R(t) ‫ܭ‬ ൌ 0

, equation S1-44 can be rewritten as

(S1-45)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. 

Using equation S1-45 and equations S1-35 or S1-36, we can write

.

(S1-46)

Both forms of equation S1-46 can be solved for N(t), but before we do this, we will elucidate the physical meaning of ‫ܭ‬ ଵ to improve our understanding of the solution.

We begin by defining a new parameter, . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint

.

(S1-48)

Using the definitions immediately above, we can write the following expression using equation S1-46:

(S1-49) Equation S1-49 can then be written as the following difference equation:

(S1-50) Equation S1-50 can also be rewritten in the following manner:

In equation S1-51, the term is independent of the infectable contacts, the change in . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. 

(S1-53) Equation S1-53 can be used to develop a recurrence relationship to explain the nature of ‫ܭ‬ ଵ .

At t = 0,

and equation S1-53 becomes

(S1-54)

We should also note that

(S1-55)

In the next time step, applying equation S1-52 to equation S1-54, we find the following expression:

Equation S1-56 can be simplified to

(S1-57)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint Repetitive application of the same logic leads to the recurrence relationship

(S1-58)

Rearranging terms and defining ݊ ∆ ‫ݐ‬ ൌ ‫ݐ‬ , we can write the following:

(S1-59)

If we now allow ∆ ‫ݐ‬ ՜ 0

, and therefore,

, and note that

we obtain the following expression for

or, more generally,

. (S1-60)

Comparing equation S1-60 to equation S1-46, we can see that . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. proportional to the strength of social distancing interventions implemented during an epidemic.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. for the remainder of this portion of the analysis. We can now find the solution to equations S1-26-28 and S1-5-7.

Using the expression for

from equation S1-46, we arrive at the following:

.

(S1-64)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint

Integrating S1-64 produces the expression for the total number of infections,

(S1-65)

Taking the limit of equation S1-65 as

, we obtain the expression for the total number of individuals who will be infected within the entire epidemic:

(S1-66)

We can write expressions for I(t) and R(t) as

(S1-69)

Equations S1-67-69 can each be differentiated and rearranged to produce these differential equations:

, and (S1-71)

(S1-72)

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) . If we also equate equations S1-71 and S1-12, we can find an expression for ߩ ሺ ‫ݐ‬ ሻ :

(S1-73)

Using these definitions for

, we can see that equations S1-70-72 are equivalent to equations S1-26-28 and equations S1-5-7. Therefore, equations S1-67-69 are solutions to the time-varying equations in the CSIR model proposed nearly 100 years ago by Kermack and McKendrick.

Equations S1-66 and S1-69 introduce an important and rather abstract concept, not previously discussed. Recall that our starting point was a reset of perspective to allow the size of the susceptible population to be variable. The susceptible population is the subpopulation of the total population that will become infected. Equation S1-69 describes precisely how the susceptible population changes as the behavior of the population, doubly exponential effect on the final number of infections that will come to be. This underscores that small changes in population behavior dramatically affect the epidemic's . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint outcome. However, the reciprocity means that the eventual number of cases produced by the epidemic is also not foreordained but rather a strong function of interventions introduced.

We gain additional insight into the meaning of the CSIR model solution by looking at equation S1-70 in more detail. Equation S1-45 can be used to rewrite equation S1-70 as

(S1-74)

The left-hand side of equation S1-74 is the rate of change in the number of new infections per person currently infected. The first two terms on the furthest right-hand side of equation S1-74, One more insight into the solution is gained from the following relationship, derived from equations S1-64-66:

(S1-75)

In words, the form of equation S1-75 is

The quantitative mathematical relationships derived from the CSIR solution in Supplement 1 characterize the dynamics of an epidemic and illustrate that strong and early intervention is critical. Equation S1-66 quantitates that the ultimate number of individuals infected in an epidemic, N ∞ , will be exponentially dependent on the number of people with which each person interacts.

The country data provide vivid examples. Both South Korea and New Zealand enacted strong and early interventions compared to other countries (3, 4) . This is reflected in their K 1 values (Table 1) Zealand was 90-99% lower than in other countries, a compelling validation of the explicit statement in the CSIR solution that strong intervention leads to exponentially more favorable outcomes.

In the USA, intervention began to have an effect around March 23 ( Figure 4B ); the number of active cases on that date was 46,136 (Table 1 ). Using the values of K 1 and K 2 in Table 1 interventions occurred widely in the USA in mid-April, well before the official reopening of the economy (13). This caused a second surge in new cases in late April and is the reason for the divergence between the observed data and the model prediction in Figure 3B . A third surge followed additional reopening activities in early fall.

As shown in Supplement 1, the CSIR solution provides an estimate of the time to the peak of new cases, t max . Using equation S1-85 and the values of K 1 and K 2 from Table 1 , the predicted peak in new cases in the USA would have occurred near March 24 if the intervention had begun on March 10. Instead, delaying effective intervention in the USA for 16 days shifted the initial peak to April 11, 16 days later, as projected; and that peak was much higher ( Figure 4B ).

As shown, too, in Supplement 1, epidemic acceleration, the instantaneous potential to change the pace of the epidemic, can be determined at any point in the epidemic and depends on the social containment actions in effect at that time (equation S1-86). What is, perhaps, less apparent, but predicted by the model, is that two countries with identical numbers of cases on a given day can, in fact, have different accelerations on the same day, and exhibit different dynamics immediately after that day.

South Korea and New Zealand (Figure 2A Table 1 ). However, since South Korea has a . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint much higher population density than New Zealand (Table 2) , it had a much higher number of interactions when the interventions were imposed and, therefore, a higher rate of acceleration as evidenced by its higher RCO at the time of intervention. Indeed, the rate of change of new cases was higher in South Korea than in New Zealand, and the later number of cases in South Korea was higher than in New Zealand ( Figure 4A and F).

Equation S1-83 clearly illustrates these lessons. As social distancing is strengthened (lower ܲ and therefore higher ‫ܭ‬ ଵ ), the replication number decreases and the epidemic slows. Early and strong interventions, especially in countries with indigenously high levels of social interaction, are necessary to stop an epidemic in the initial stages; and reopening actions, enacted too early, can reignite the epidemic, dramatically increasing the number of cases. The astonishing magnitude of the effects, driven by only a few days of delay, derive from the doubly exponential nature of the underlying relationships.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint SUPPLEMENT 3

Ending an ongoing epidemic

The CSIR solution can also be used to design measures to end an epidemic in an advanced stage.

The management plan is built by first using equation S1-90 to predict how many days a given level of intervention, ‫ܭ‬ ଵ , takes to reduce the new daily cases by a target fraction:

is the time to the desired reduction and ‫ܦ‬ ௧ is the target as a fraction of the current level of new cases per day.

For example, if a country were to target a 90% reduction of new cases per day (e.g., from 50,000 (7)).

Since భ ൌ ܲ , K 1 = 0.2 characterizes a lockdown in which people in the country can each only have one plausibly infectious contact with a little over one specific person for the containment duration. This does not mean they cannot contact anyone other than the one person; but they . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint must use care, using masks and proper distancing, to ensure there is no plausibly infectious contact with anyone other than the one person.

Returning to the planning example: once the initial 90% reduction is achieved, a reasonable next step might be to relax the social containment to a level that allows the economy to remain viable while preventing the epidemic from erupting again. The level of ‫ܭ‬ ଵ necessary to achieve a chosen target can be again found using equation S1-90. If an additional 90% reduction in new cases per day is desired, and a period of 90 days is tolerable for that reduction, then the new level of ‫ܭ‬ ଵ needed is ~0.025. This equates to allowing each person to be in contact with 7 specific people, in an infectable way, for 90 days. Note that this is 3 times less stringent than the original US shutdown level in April 2020. Thus, with a well-planned approach, a country can reduce its new daily cases by 99% in approximately 100 days, enabling the country to control, and essentially end, the epidemic while maintaining economic viability. , low enough to prevent an outbreak and yet continue decreasing the new cases per day, can be found using equation S1-87.

The progress of interventions is easily monitored using the RCO, as the curve for South Korea illustrates ( Figure 2A ). Had this country maintained the implemented level of distancing measures, the data would have followed the initial slope. However, the actual data departed from the slope, heralding the failures in (or relaxations of) social distancing, which were later documented to have occurred during the indicated time frame (3) (circled data, Figure 2A ).

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint

A contribution to the mathematical theory of epidemics

Early draconian social distancing may be suboptimal for fighting the COVID-19 epidemic

South Korea's health minister on how his country is beating coronavirus without a lockdown

New Zealand isn't just flattening the curve. It's squashing it. The Washington Post

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. Equations S1-76 and S1-77 illustrate the logic of the CSIR solution in terms of probabilities. (S1-78)

The solution to equations S1-5-9 was developed assuming that However, both assumptions can be tested independently from the country correlations.The check of the assumption that ‫ܭ‬ ் is a constant can be made by substituting equation S1-63into equation S1-65 and solving for ‫ܭ‬ ் t. Doing this, we find the following expression:. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) , then we can also write this expression:is a constant, then equation S1-80 predicts thatis a linear function of time.Excepting ‫ܣ‬ ଵ , all the quantities on the left-hand side of equation S1-70 can be derived from the country data. Therefore, we can estimate the value of ‫ܭ‬ ் using equation S1-79 (or S1-80) to find the value of ‫ܣ‬ ଵ that best fits a straight line. This fit is illustrated in Figure 5 and shows a correlation coefficient of 0.955, which strongly supports the assumption that ‫ܭ‬ ் is a constant.Equation S1-64 can be rewritten to define an important relationship:

is the time between the start of the epidemic and the time at which the first data point in the epidemic is measured.The RCO expression is convenient in that it transforms equation S1-64 into an equation linear in ‫ܭ‬ ଵ . We call this expression the "rate of change operator" since it is the rate of change of new cases scaled by the current daily case number. Equation S1-81 is the equation that was fit to the . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint country data in Figure 2 , and used to find the parameters in Table 1 . Note that if ‫ܭ‬ ଵ decreases and if it approaches zero, the growth of the epidemic will be nearly exponential.Similarly, we can use the fact that equation S1-46 predicts that the RCO measure is proportionalwhere ‫ܭ‬ ଵ is a measure of the interactive behavior of the population. If an independent measure of people's mobility during the epidemic should also be linearly related to the RCO, we can have additional confidence in the veracity of the CSIR solution. Google has compiled different measures of people's mobility derived from mobile phone data (11). One of these measures is termed the residential mobility measure (RMM). The RMM is a measure of the percent increase or decrease that people stayed in their residence during the pandemic relative to a baseline measured over 5 weeks starting on January 3, 2020. Since To test this, we plotted the daily integral of the RMM for the six countries we analyzed, against the daily RCO. These plots appear in Figure 6 and the hypothesized linear relationship is clear.

Substituting the values for R(t) and N(t) from equations S1-68 and S1-65 into equation S1-73, we can arrive at an expression for ߩ ሺ ‫ݐ‬ ሻ as a function of time:(S1-82). CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint From equation S1-27, we can see that the number of infections, I(t), will begin to decrease when , we can write the following criteria for when the epidemic will begin to decline:(S1-83)Using equation S1-83, we obtain the following expression for when the decline begins:(S1-84)If we differentiate both sides of equation S1-64, we obtain an expression identical to equation S1-84: Another important expression is the rate of acceleration of the epidemic:. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint. (S1-

Equation S1-86, with its four equivalent expressions, is a demonstration of the power that an . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) (S1-88) and using equations S1-64 and S1-65, we arrive at the following expression:and we can derive equation S1-90 from the remaining terms:(S1-90). CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. , that a level of social containment, ‫ܭ‬ ଵ , will be required to achieve a fraction of daily cases, ‫ܦ‬ ௧ , compared to the current level.. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. and equation S1-65 asDistancing measures tend to be constant for many days at a time, so for this analysis, we assume . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted March 5, 2021. Equations S3-1 and S3-2 can also be rewritten as, and (S3-8). CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. Using equation S1-64, we arrive at the following expressions:( . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) If new infections are introduced into a portion of the population that has been thus far disconnected from the previously infected area, then the assumption of contiguity has been violated. This is a common situation when infected people travel from an infected area into an area that was previously uninfected or had not yet seen significant numbers of infections. We label this a Type 2 outbreak. Equation S1-65 must be modified to predict the number of cases in an epidemic that is having , where ‫ݔ‬ denotes the number of the outbreak, equation S1-65 can be written as:If a new outbreak occurs in a previously unaffected area of a country, then equation S1-65 can be modified as follows: . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)The copyright holder for this preprint this version posted March 5, 2021. , and ܰ ௫ need to be determined independently.While an epidemic is underway, a Type 2 outbreak can be detected by monitoring the slope of the RCO curve. If a positive slope is detected in an RCO curve, a Type 2 outbreak has occurred. This is an indication that immediate action, within days, is required from policy makers to strengthen intervention measures and prevent the outbreak from overwhelming prior progress in controlling the epidemic.If the disease changes its transmissibility through mutation, this can also be detected by monitoring the RCO curve. In this situation, a proper fit of the parameters in equation S1-71 will not be possible and a modification of ‫ܭ‬ ் will be required to accommodate the change.. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)The copyright holder for this preprint this version posted March 5, 2021. ; https://doi.org/10.1101/2021.03.02.21252781 doi: medRxiv preprint