key: cord-0980437-6l56kdr4
authors: Cherednik, I.
title: Modeling the waves of Covid-19
date: 2021-06-22
journal: nan
DOI: 10.1101/2021.06.15.21258969
sha: 343d295c215bf3b3ada3e9346dfb142b8d4fa7d6
doc_id: 980437
cord_uid: 6l56kdr4

Classical approaches to modeling the spread of epidemics are based on two assumptions: the exponential growth of the total number of infections and the saturation due to the herd immunity. With Covid-19, the growth is essentially power-type, especially during the middle stages, and the saturation is currently mostly due to the protective measures. Focusing on these features and the role of epidemic management, we obtain differential equations for the total number of detected cases of Covid-19, which describe the actual curves in many countries almost with the accuracy of physics laws. The two-phase solution we propose worked very well almost for the whole periods of the spread practically in all countries we analyzed that reached the saturation during the first waves. Bessel functions play the key role in our approach. Due to a very small number of parameters, namely, the initial transmission rate and the intensity of the hard and soft measures, we obtain a convincing explanation of the surprising uniformity of the curves of the total numbers of detected infections in many different areas. This theory can serve as a tool for forecasting the epidemic spread and evaluating the efficiency of the protective measures, which is very much needed for epidemics. For instance, the accuracy was very high for the early projection for the 3rd wave in the USA; projections for India, South Africa and UK are considered.

which approach results in a very exact modeling of Covid-19. This is actually similar to [Has] : "The question of interest was the time course of the epidemic, rather than the final state, which is always one where the disease dies out". The "predator-pray" system for us is when the protective measures (including self-imposed ones) play the role of "predator", and the "pray" is the number of infections. This is different from the SIR-type models, applicable mostly to the initial periods of exponential spread and to final stages of epidemics. The SIR model was suggested in early 20th century. Since then, it was developed, but the exponential growth until the herd immunity is approached remains its essential feature. As we will show, the asymptotic periodicity of Bessel functions is absolutely relevant here. Generally, Bessel function and processes are important in mathematics and physics, but they were not employed for epidemics and in invasion ecology as far as we know.

Our approach seems promising in ecology. More specifically, Bessel functions can presumably describe various continuous 2-species models. Following [Has] , the discretization, different time-scales, 3-species models are natural further steps. See also [HL, LPP] . Generally, basic hypergeometric functions and their variants are expected to occur.

The discretization will be discussed only a little. We also do not consider in this paper the concept of Momentum Risk Taking from [Ch2] , somewhat similar to Kahneman's "thinking-fast", which can be considered as a behavioral counterpart of the "transiencies" in ecology.

The most ambitious here are the expectations that the same ODE model the processes of momentum decision making in our brain, but this is very preliminary. The number of neurons involved in the "momentum" analysis of some event is restricted here by the "predator", the expected allocation of (very limited) resources of our brain for this particular task. The asymptotic periodicity of Bessel functions sets here some limits. Generally, it is surprising that Bessel functions, invented by Daniel Bernoulli long ago, are not one of the main tools in mathematical theory of epidemics, ecology, behavioral science, and beyond. Hopefully this will change.

The usage of the basic and current reproduction numbers R 0 , R is common for epidemics. The basic one, R 0 , is defined as the initial average number of people infected by one person who contracted the virus; see [CJLP, Co, CD, DHB, He, HL] . However, R can be used only qualitatively for Covid-19 and other epidemics of power growth: the formula u(t) ∼ const R t for the total number of infections will stop working very quickly and cannot be of real help for the forecasting the spread of epidemics without very significant corrections.

Even R close to 1 would quickly begin to contradict the actual growth of infections of Covid-19. However, it is not unusual when they were reported as 2, 0.7, or so in the middle stages of the waves of Covid-19; they are provided constantly by Robert Kox Institute (for Germany) and other centers. This is questionable to us.

One of the possibilities to adjust SIR to the power growth of Covid-19 and other epidemics of non-exponential type (there were many such), is to assume that R ∼ 1 and that it is non-dominant, i.e. to invoke the theory of resonances. This provides a polynomial growth of the total number of infections, but such models are unstable mathematically. Our modeling is different: a combination of the "local herd immunity" with the role of active management.

One of the most efficient protective measures is (and always was) self-restriction of our contacts. Here the disease control centers are supposed to provide current information on the spread of the epidemic. This is of course combined with quite a spectrum of other measures.

It was already intensively discussed in the literature that the herd immunity can influence the spread of Covid-19 well before it reaches the levels of 70% or so. See e.g. [BBT] . The protective measure play a significant role in this reduction. However their relaxation can result in the recurrence of the waves of infection.

We actually make the next step: claim that local herd immunity begins almost from day one and that it results in the power-type growth of the total number of infections: u(t) ∼ t c . Here the exponent c depends very much on the time from the beginning of the wave of infections, which is addressed in our approach via Bessel functions.

Our theory was posted in the middle of April, when the saturation of the spread was observed only in several countries; they were mostly in phase one, under mode (A) in our terminology. We also provided a variant for phase two, under mode (B), when the hard measures are significantly reduced. The (B)-mode system of ODE appeared really necessary for correct modeling the spread. Qualitatively, phase two is the switch to some less aggressive management due to relatively low numbers of daily infections, including restricting and self-restricting our contacts.

For the initial growth u(t) ∼∼ t c of the number of detected cases, the second phase is described well by u(t) ∼∼ t c/2 cos(d log(t)) for some d. The passage from the Bessel-type functions u(t) for phase 1 to such ones can be clearly seen in many countries. Though the Bessel-type u(t) worked well almost till the saturation in quite a few, including our forecast for the 3rd wave in the USA.

The spread of Covid-19 in the USA was mathematically quite a challenge for us; the results of our prior efforts are systematically reported in [Ch1] . The first wave in the USA went through several stages, more than with any other countries we considered. Our understanding is that it was so mostly because the hard protective measures were constantly relaxed in the USA on the first signs of improvements, well before the actual saturation. This is in contrast to Europe and several countries in Asia. It was somewhat similar in UK, but it eventually reached phase 2 and the saturation of its first wave.

The costs and consequences of hard measures, especially lockdowns, are huge. However the saturation due to the hard measures is of unstable nature and the recurrence of the epidemic is quite likely if they are reduced or abandoned. Our theory generally provides the way to evaluate the efficiency of protective measures and employ them properly, but this is quite a challenge even if advanced mathematical means are used. See e.g. [FRAF] .

Prior approaches. There is increasing number of works where the power growth of the total number of infections is considered for modeling Covid-19. Let us mention at least [Ch1, MBS, MH, TKH] .

Let us begin with [CD] (well before . An ambiguity with the definition and practical calculation of R 0 , R is discussed there: "It is reassuring to know, however, that the sign of R 0 −1 is independent of the decomposition used and that the prediction of exponential growth or decay is therefore correctly made by any of the counting schemes." This is our impression too: the sign of R−1 is what is mainly used practically, not the exact value of R (calculated by some simple formulas). We note that the spread is mostly assumed of exponential type in [CD] . Let us quote: "As far as we know, little can be said in general about the exceptional case that R 0 is not strictly dominant".

In [MBS] , the authors comment on the power growth of the spread of Covid-19: "the nature is full of surprises". In [TKH] : "this new contamination regime is hard to explain by traditional models". In our one: "power law of epidemics must be the starting point of any analysis if we want our mathematical models to be up to date". See also [Ray] and works mentioned there concerning a potential usage of small-world interaction network, where individuals are assumed to contact (mostly) local neighbors and have occasional long-range connections.

In a different direction, paper [BBT] and some other works suggest that the levels of herd immunity sufficient to impact the spread of Covid-19 can be significantly lower than the "classical" 70% or so: as low as 40% in some areas due to the population heterogeneity. From this viewpoint, we make the next step in this direction, which seems quite natural. Our starting assumption is that local herd immunity shapes the spread from the very beginning of epidemics and quickly reduces its exponential growth to the power one. This is related to the concept of small-world. Spatial modeling. This approach is actually similar to the one via "small-world". The graph of contacts, especially geographically related ones, is the key for spatial modeling. For instance, Fig. 8 in [BP] shows the initial spread of Covid-19 in Germany, which is naturally related to the geographic locations. The authors change there the SEIRD model (susceptible-exposed-infected-recovered-deceased) to SEIQRD by adding the quarantined compartment. Obviously quarantines and travel restrictions are important to model epidemics. See also Fig. 2 from [KBLK] (Germany too). The population heterogeneity is an important consideration here.

The challenge is to produce meso-scale forecasts, with specific information and for some concrete locations. See e.g. [JCO] (for Texas). The total number of detected cases and other general data for the whole country are insufficient for those in charge of practical managing the epidemic. Partial differential reaction-diffusion equations are used in [JCO] ; see also [VLAB] .

Classical methods "least-squares", "Bayes", "k nearest neighbors", and various statistical tools can be used for "meso-forecasts". Deep machine learning is a possibility here (which we will not touch upon). AI-systems generally do not help much with theoretical understanding the processes. There are too many parameters and the uniqueness of the optimal ones is not granted. However they "almost always" provide satisfactory projections. ODE, PDE, SDE are not strictly necessary for machine forecasting, but they are generally in terms of meaningful parameters, which can be of importance theoretically.

Surprisingly, the dimension of the graph of contacts is basically sufficient for forecasting the number of infections in sufficiently large areas: states, regions. This is our parameter c; the total number of infections is then ∼ t c at some early stages. This may hold for small areas, and even for the viral load in infected individuals (at the middle stages), but this is a macro-phenomenon in this paper. The uniformity of the curves of total numbers of detected infections in so many so different countries is still mysterious to us, though we think that we found proper mathematical tools to address this.

Power growth. The main problem with modeling Covid-19 appeared actually not the power growth itself, "power law of epidemics", a starting point for us. This alone is insufficient for forecasting. Understanding the saturation is the key challenge for modeling epidemics. However, we must note here that the power growth of the number of infections is not commonly accepted, in spite of ample evidence during the Covid-19 pandemic.

Practically, the exponent c and the corresponding scaling coefficient C heavily depend on the time passed from the beginning of the corresponding wave of the infection. An exact mathematical model of this time-dependence is necessary here, which was proposed in [Ch1] using the Bessel functions.

For instance, the approach of [MBS] to the power growth was of experimental kind. The corresponding exponents depend very much on the considered periods. So the data in Figure 1 in [MBS] and in similar papers mainly show qualitatively that the growth is no greater than power. However the latter is quite convincing! The parameter c we operate with is quite different from their exponents; it is for the whole period of the wave, even including the 2nd phase (where no Bessel functions are used). We mention that mathematically, the authors used the SEIR model (susceptible-exposedinfectious-removed), which does not result in the power growth; however, "small world" is mentioned there as a possibility.

Paper [TKH] is based on the Poissonian small-world network. This approach results in the linear growth (c ≃ 1) of the number of infections. The linear growth is indeed present indeed near the turning points of the curves of total numbers of infections, but it is far from . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 linear anywhere apart from these middle portions. The explanation of the linear growth and the saturation in [TKH] is very different from what we proposed in [Ch1] and present in this paper.

Among many confirmations of the power growth of the total number of infections of Covid-19, the period 3/20-10/7 (2021) in India is very convincing; see Fig. 9 . Here u(t) = const t c for c = 3.65 is practically exact (almost without modifications) for the total number of detected cases in India for a very long period: for at least 5 months (!).

In this figure, the main parameters were determined on 08/03. This forecast was posted on 10/07 (2021); to expire on 11/06 (the maximum of the u-function). It matched very well the actual curve of detected cases. This was based on our usage of Bessel functions. Here and almost in any countries, some linear-type growth is expected after the top of the Bessel-type curve u(t), which can be seen in the graph for India after 11/06. In our theory, this period is described by the formula for mode (B); this is phases 2 in our terminology (for any waves).

Saturation due to hard measures. For us, the saturation, followed by some linear-type growth of the total number of infection, is due to protective measures, mostly the hard ones. They are those imposed by authorities in charge, but self-restrictions are equally important here.

The key is detection-isolation-tracing, which includes closing the places where the spread of infection is the most likely. The societal cost of hard measures is huge, but they proved to reduced the spread efficiently. Mathematically, the vaccination is of this type too, a very hard measure by any standards.

It is not disputed that the saturation of the first waves of Covid-19 was not due to the herd immunity. The latter probably requires about 40%-60% (smaller than "classical" 70%) of all susceptible population to be infected and recovered [BBT] . It was far from these levels during the first waves. Thus, the saturation mechanisms of SIR-type models are not applicable here, at least for the 1st and 2nd waves.

The timing and the intensity of the second waves clearly confirm the validity of our approach to modeling, based on the prime role of protective measures, mostly the hard ones. Recurrence of epidemics is quite frequent; see e.g. [HL] . However the second waves of Covid-19 began unusually quickly (for epidemics), sometimes even on the top of the unfinished first waves, as it happened in the USA. The relaxation . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 of hard measures closer to the end of the first and further waves seems the only logical explanation of such sequences of waves.

In 2020, the summer vacations (and closed schools) in Western Europe were a clear instance of (broadly understood) hard measures. They did reduce the spread. However, at the end of August, the second waves began practically everywhere in Europe. Similarly, the number of new detected infections began to grow quickly (again) in the USA from the middle of September (2020).

Mathematically, our exponent "c", which we call the initial transmission rate, appeared increasing from the first to the second waves in many countries. This parameter is one of the main on our theory; it combines the virus transmissibility strength and the number of contacts in the areas at the beginning of the current wave of the infection. Thus, by reducing the protective measures, c can be expected to be essentially back to its levels in the beginning of the epidemic, though with a tendency to increase. The increase of c from the 1st wave to the 2nd and the 3rd is related to the general reduction of Covid-related restrictions (including our own behavior). The new strains of Covid-19 contributed too; its evolution is of obvious importance.

The second key parameter of our theory, the intensity a of protective measures, dropped very significantly for the second waves, which was expected. The intensity of the hard measures during the 1st wave was difficult to sustain during the next waves.

Qualitatively the duration of the wave is 1 √ a ; quantitatively, Bessel functions must be used here for exact modeling. So, mathematically, we essentially repeat the first waves, but now with significantly lower levels of hard measures, longer periods of intensive infections, and higher magnitudes of the curves of total numbers of detected infections. Though, the example of 2 waves in India is somewhat exceptional.

Power Law of Epidemics. With such complex processes as epidemics, there can be of course multiple factors contributing to the power growth, biological ones included [CLL] . The "justification" from [Ch1] goes as follows. First, we assume that infected people mostly transmit the disease to their (susceptible) neighbors, and that the population is distributed uniformly. The second assumption is that the wave of the infections expands linearly in a proper graph of contacts.

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 The third assumption, the principle of local herd immunity, is that people "inside the infection zone" do not transmit the disease because they are surrounded by those already infected or recovered, i.e. the border of this zone mostly contributes to the spread of this disease. This readily gives that u(t) ∼ t 2 or greater (in the absence of protective measures). Indeed the lowest c we observed was c = 2.2 (the 1st wave in the USA).

People from infected zones do shopping, travel, visit friends. So the higher dimensions are needed to imbed the graph of contacts into some R N providing that the geometric distances between points representing people are essentially the numbers of links between them, i.e. that these distances reflect the intensity of the contacts.

Upon this embedding, we assume the uniform distribution of the points in R N representing people, and the linear spread of the disease in R N . This is basic physics. Then, indeed, u(t) ∼ Ct c , where c is the "dimension" of the image of this graph, a number from 2 to N .

Next, we represent such u(t) as a solution of the differential equation du(t)/dt = cu(t)/t. This is standard when we need to add "external forces", which will be protective measures. We argued above that the exponential growth is generally unsustainable. However the power growth is unsustainable long term too. This will be "corrected" as follows.

Adding protection. Combining the initial power growth of the total number of detected infections u(t) with the impact of protective measures we obtain the following two systems of differential equations:

Here t is the time from the beginning of the intensive growth of infections, not always the very beginning of the corresponding wave of Covid-19 but sufficiently close to it. System (1) describes the impact of hard measures under the most aggressive response to the spread. The second system describes the impact of the soft measures: some travel restrictions, wearing the protective masks and social distancing are typical. We called these two modes (A) and (B) in [Ch1, Ch3] .

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101/2021.06.15.21258969 doi: medRxiv preprint When a = 0, d = 0, we obtain the power growth u(t) ∼ Ct c ; so c can be measured experimentally during the initial stages of Covid-19 and is supposed to be the same for (1) and (2). Mostly it was in the range 2.2 ≤ c ≤ 2.8 (wave 1), but reached c = 4.5, 5.5 in Brazil and India.

There is a variant of these systems, when the second equation in (1) is replaced by that from (2), called the transitional (AB)-mode in [Ch1] . It modeled reasonably the spread in the USA, UK, and Brazil, but the usage of (1) and (2) in our two-phase solution appeared sufficient for many countries without mode (AB).

The protection function p(t) for (1) is basically the number of prevented infections. More exactly,

where the sum is over all infected individual isolated at the moments 0 < t i < t for some constant a, the intensity of "isolations". We assume that if not isolated, this group of people would contribute p(t) to du(t)/dt (the transmission rate is taken 1 for them). For (2), p(t) ∼ (the number of infected people wearing the masks before t), and, similarly, for other (self-)restrictions.

SIR-type modification Let us touch upon the modification of system (2) under the assumption that u(t) is bounded. Assuming that u(t) < 1, we multiply the right-hand side of the 1st equation by (1 − u(t) ), which models the interaction of infected individual with the remaining (susceptible) ones. We obtain:

In the absence of p(t), it is a well-known logistic equation, with the following modification: the interaction is proportional here to 1/t. The 2nd equation remains unchanged. The resulting system can be solved numerically. It is not clear whether the corresponding solutions are more relevant than those for the original system (2). For the following modification of the 2nd equation, this system can be readily integrated: dp(t)/dt = b d u(t)/dt. One has:

Related processes. Both systems are actually from [Ch2] , where they were used to describe the dynamic of the (relative) stock prices . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101/2021.06.15.21258969 doi: medRxiv preprint p(t) under news driven momentum trading. The function u(t) there was the news propagation triggered by some event. It is of power growth in terms of time t passed from the event, but the exponent c is generally significantly smaller than 1, especially for short-term trading.

The arguments there were from behavioral finance. This is actually related; the behavioral aspects of epidemics are of obvious importance [St] . However financial news fades, and this happens quickly; this is very different for the spread of epidemics. System (1) described in [Ch2] profit taking in stock markets; the second one modeled the "usual" news-driven investing.

As a matter of fact, these two systems are of very general nature. For instance, they are supposed to occur in any momentum risk taking. This concept, MRT for short, is from [Ch2] ; it is somewhat similar to Kahneman's "thinking-fast" [Ka] . Managing epidemics on the basis of the current data is very much momentum. As in stock markets, people and authorities in charge must react promptly to any change of the situation. Another example there is tree growth, though there were no p(t) and the arguments were somewhat different.

It was expected in [Ch2, Ch1] , though without biological evidence, that both systems of equations may describe real neural processes in our brain. Here u(t) is the number of neurons involved in the analysis of a particular event at the moment t, counted from the event, and p(t) is the expected importance of this event vs. other ones and the corresponding expected brain resources needed for its analysis. I.e. p(t) is basically the expected allocation of resources, which are very limited in our brain. We do not know much about the ways our brain work, but the confirmation of the power laws and related saturations are solid in the stock markets and, as we demonstrate, in epidemics.

We note that a significant part of [Ch2] is devoted to the discretization. Decision-making always requires some action potentials, i.e. it is discrete by its nature. With epidemics, the usage of ODE worked very well so far, though potentially the discretization can become important for our approach too.

Toward discretization. We begin with some basics. Let u n be the total number of infections at the nth moment from the beginning of some epidemic. Infected individuals transmit the disease mostly during some initial period, which we will make the unit of time. Let it be 1 week . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10. 1101 (indeed, about 1-2 weeks for Covid-19). The recurrence relation and the corresponding quadratic equation are:

Here R 0 is the initial reproduction number.

We obtain:

If C 2 ̸ = 0 and R 0 > 1, the growth of the total number of infections will be exponential. This is unless the herd immunity is expected to be reached, when a difference version of SIR is needed. We will not address this here.

Power growth. Basically, there are 3 main mathematical possibilities to ensure a power growth:

(i) the presence of "predators", forces restricting the epidemic spread, for example, various protective measures, (ii) when R approaches 1, which results in "resonances" and can potentially provide some linear growth, (iii) when the "birth rate", the transmission rate in this context, becomes inversely proportional to time.

Obviously (i) is applicable: we do fight epidemics. The resonances and linear growth occur when R ≈ 1, but this is unstable and does not seem of actual importance for modeling epidemics. We think, a combination of (i) and (iii) is the key in epidemics.

The spread of any disease is the growth of the "circles" of those infected; these are combinatorial circles, not geometric ones (in a map of the affected area). The contacts of infected people are not only with their immediate neighbors; people work, study, do shopping, travel. This is the concept of "small world" and the basis of spatial modeling. Our main assumption is that the rate of change of these circles can be expected constant. However they are not planar ones.

The individuals at the frontier of such a circle transmit the disease the most because:

(a) they are the "latest" and therefore in their most infectious stages, (b) people inside the circle are "surrounded" by those with immunity.

We assume that the infected people "inside the circle" contribute to the transmission of the disease significantly smaller than those at its boundary. Then the circle of infected people can be presented as a ball of dimension c; accordingly, the growth of u n will be ∼ n c . Here c can . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101/2021.06.15.21258969 doi: medRxiv preprint be any positive number, not only an integer. For instance, c ≈ 2 if our contacts are mostly with those who live near us, but c ≈ 3.65 in Fig. 9 . We never saw countries with c < 2 during the 1st and the 2nd waves. Recall that c, the initial transmission rate, reflects the virus transmissibility strength and the number of contacts in the area at the beginning of the wave.

Disregarding "maturity", a period when people are already infected but not infectious, the number of "newly infected people" is basically the area of the boundary of this "ball", i.e. the area of the corresponding sphere. Presenting this number as c u n−1 /(n − 1), we arrive at the following recurrence: u n − u n−1 = c u n−1 (n−1) . Asymptotically, n c is a solution of this recurrence. Let us comment on this.

We have

Tree growth. We mention that the equation u n = u n−1 + u n−2 (n−2) describes middle stages of tree growth; here the time-unit is naturally 1 year and the "maturity" is set to 1. Its obvious solution is u n = n. Under the initial conditions u 1 = 0, u 2 = 1, the corresponding solution tends to n/e, directly related to derangements D n in combinatorics: u n = D n /(n − 1)!. We omit the experimental support for this model.

The following recurrence is convenient to model the saturation stage of tree growth: u n = u n−1 + (n−1) (n−2)n(n+1) u n−2 . Its solution with the initial conditions u 1 = 0, u 2 = 1 3 tends to 2J 0 (2) ≈ 0.447782 as n → ∞, for J 0 (2) = 1 − 1 + 1 (2!) 2 − 1 (3!) 2 + . . . ; u n = n n+1 is another solution. Adding protective measures. In the most aggressive variant, we have:

where a is the intensity of "hard measures", p n is essentially the number of people protected from the virus due to the measures until n.

Basically, we subtract here the number of "preventions" from the total number of infections. In the second relation, the protection is assumed proportional to the current number of cases. This is parallel to the differential case. The main point here is that the "isolation"

. CC-BY-NC-ND 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) Next, we eliminate p n from these relations. We obtain that p n−1 = (1 + c n−1 )u n−1 − u n , p n = (1 + c n )u n − u n+1 , and

Finally, we obtain the recurrence relation:

When a = 0, we have the solution u n = n for any c. This is for u 1 = 1, u 2 = 2. For a = 0 and arbitrary initial conditions u 1 , u 2 , the "function" u n is essentially proportional to n c , which follows from (3). When c = 0, the equation is not applicable: p n = u n − u n+1 , the protective measures, cannot be negative. Fig. 1 is an example of the calculation with (4), where a = 0.1, c = 2.2, and we begin with u 1 = 1, u 2 = 3. We must stop at n = 8 (week 8), which is the saturation; the total number of cases cannot decrease.

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10. 1101 Near and after the maximum of u, a different recurrence must be used, which we omit.

Two-phase solution. The solutions of (1) and (2) we need are

Here

it will be used to model later stages of the waves of Covid-19.

Our two-phase solution is the usage of a proper linear combination of u 1,2 for phase 1, till its saturation, and then the usage of u B for phase 2. It proved to be quite exact for modeling the curves of total numbers of detected infections of Covid-19. For t ≈ 0: u 1 (t) ≈ t c and u 2 (t) is approximately ∼ t. I.e. u 1 dominates; it is the key for forecasting.

The second fundamental solution of system 2 is with sin instead of cos. We note that when the protective measures are modest, we obtain

i.e. t c in the beginning of the spread diminishes to t c/2 and then remains unchanged. This is of importance, but we will not touch the range D > 0 in this work.

The following examples of the 1st waves in 2020. mainly follow [Ch1] . Italy: 2/22-5/22. Figure 2 . The starting point was 2/22/2020, when the total number of infections was 17. We subtract this initial value here and below when calculating our dots, the total numbers of detected infections. One has: (1, t) )), c = 2.6, a = 0.2, d = 0.5.

We use here both fundamental solutions u 1,2 (t) of system (1).

Germany: 3/07-5/22. See Figure 3 . We began with the initial number of total infections 684 (subtracted). This was approximately the moment when a systematic management began. One has: (1, t) )), a = 0.35, d = 0.56.

. CC-BY-NC-ND 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-ND 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. Japan: 3/20-5/22. See Figure 4 . There was some prior stage; we subtract 950, the total number of infections on March 20. The curve for Japan is not too smooth, which is not unusual. However it is managed well by our 2-phase solution : (1, t) )), a = 0.3, d = 0.6.

UK: 03/16-06/13. This country was a challenge for us, though it "eventually" managed to reach phase 2. Actually the red dots for UK are modeled better with the transitional (AB)-mode. However, we prefer to stick to the "original" u(t) determined for the period till April 15. The two-phase solution is a combination of two phases separated by a linear period, about 10 days. See Figure 5 . The formulas are: (1, t) )), d = 0.465.

The Netherlands: 03/13-5/22. The u-function here is with the same a, c as for UK. The parameter d = 0.54 is different from that for UK . CC-BY-NC-ND 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-ND 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 June 22, 2021. ; https://doi.org/10.1101/2021.06.15.21258969 doi: medRxiv preprint (d = 0.465). This could be expected; the process toward the saturation of phase 2 was slower for UK.

See Figure 6 . The number of the total case was 383 on 3/13, the beginning of the intensive spread from our perspective. The usage of the dominant u 1 appeared sufficient: (1, t) )), d = 0.54.

The mathematical similarity of the second and 3rd waves to the first waves is very remarkable, a strong confirmation of our approach. The parameters a, c, b though change, which is generally in a quite understandable way. We will begin with the 2nd wave in the Netherlands, from 08/24/2020.

The Netherlands: the 2nd wave, 2020. The second waves were quite uniform in Western Europe. The Netherlands is convenient to demonstrate the evolution of our parameters, because the corresponding ufunction does no involve too much of the second, non-dominating, Bessel-type solution. Generally, both are present.

The similarity of Fig. 7 and Fig 6 is obvious. This is similar to the qualitative similarity of the 1st, the 2nd, and the 3rd waves in the USA to be discussed below.

The parameter c significantly increased in the Netherlands: from 2.4 (the 1st wave) to 3.4 (the 2nd). The intensity of the hard measures understandably dropped: from 0.2 to 0.085. Such changes are actually common for the second waves in Europe. The parameter d = 0.43 diminished from 0.54. Our projection worked well until the beginning of December, actually until the pause between the 2nd and the 3rd waves.

In December, almost all Western Europe switched to a linear-type growth of the total number of detected infections. Later, the 3rd waves there on top of the unfinished 2nd waves began there. The usage of protective measures became less stable, which can be partially due to the holiday season and the new strain of Covid-19 was on its way to dominance. Anyway, the intensity parameters a diminished anywhere in Europe and the USA vs. those for the 1st waves.

. CC-BY-NC-ND 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. For the second wave in the Netherlands, one has:

u B (t) = 1.65 t c/2 cos(d log (M ax(1, t) )), a = 0.085, d = 0.43.

In Europe and the USA, there were significant fluctuations of the spread during the winter. Some hard measures were re-introduced in October-January. The vaccinations began in Spring, certainly very hard measures. This made and will make us closer to the herd immunity (for these particular strains). The self-imposed restrictions are of course of obvious importance here; their reduction can be one of the reasons for the increase of the initial transmission rates c during the 2nd waves vs. those for the 1st ones.

Japan: the 4th wave (in process). It began in March, 2021. We take the period 03/15-5/25, 2021. The formula for u(t) is:

The tendency for c to somewhat increase and for a to drop is similar for that in other countries we considered. Here we did not create the control period; all dots are red. It is clear from our modeling that the . CC-BY-NC-ND 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. 1st phase of the 4th wave approaches its end in Japan, and there is a switch in process to the 2nd phase.

The 1st wave in India: 3/20-10/07-11/20, 2020. This country provide important mathematical patterns of the dynamic of the spread of Covid-19. The (clear) first wave was later than in quite a few countries, but it was with the greatest c we observed. The starting number of detected cases was 191, which was subtracted.

The power function 0.0125(t + 0.07) 3.65 gave a surprisingly good approximation for more than 5 months; see Fig. 9 . The parameter c began to decrease much faster in other countries. Obviously, the size of population is a factor here. This unusual stability the exponent can be also linked to a relatively low level of the active management in this country during the first wave. Of course the density of the population (which is huge in many parts of India) and the general number of contacts are very significant factors too. This is in India and anywhere.

Such a clear power growth is a very convincing argument in favor of the power law of epidemics. It is like this can be seen in all countries. It was mistakenly called exponential in some publications, which was of course not the case.

. CC-BY-NC-ND 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. In Fig. 9 , the parameters were determined around 08/03. Generally, the parameters determined before the turning point of the curves of total cases must be considered conditional, though we did it successfully in several cases before these points. The forecast posted on 10/07 was that the curve of the total number of detected infections would reach its technical saturation on November 6 with the number 8.25M of the cases. It matched the actual number of cases well.

As always, a linear-type growth (mode (B)) is expected after the top of the Bessel-type curve u(t), which can be seen in the graph. This period is described by u B . The parameter c is the same for phase 2 as it was for phase 1 (in the Bessel-type u(t)). Also, u B begins from the starting point of the wave: "as if it were no 1st phase".

We omit u B for India. Let us mention here that our computer forecasting programs are for the 2nd phase only. They find the best u B (from point 0) approximating the last 20 points. See below.

Here y =cases/10K; similarly, y is the total number of cases divided by proper powers of 10 in the other charts we will consider. Say, divided by 100K for the USA. The x-axis is always time in days from the beginning of the curve. The red-blue-black dots give the corresponding actual total numbers of the detected cases. The u-function for India is:

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101/2021.06.15.21258969 doi: medRxiv preprint where c = 5.75, a = 0.035. It matched very well the actual numbers of cases till the middle of December (till the 2nd phase of the 1st wave).

The 2nd wave in India: 03/25-05/25, 2021. The corresponding ufunction is

, where c = 2.85, a = 0.14.

We note that the parameter c, which is the key, became smaller than during the 1st wave in India, though there is a somewhat greater contribution of the second solution u 2 , the term 0.4J 0.5−c/2 (·). Generally, this term affects C, but does not influence much c, a. Its "role" is to adjust better the early stage of the wave. At later stages J 0.5−c/2 becomes negative. This is why C must become greater in its presence with positive coefficient; this term was 0.2J 0.5−c/2 (·) for the 1st wave.

The coefficient a for the 2nd wave is quite similar to those for the 1st waves in Europe, better than those for further waves in Europe and the USA. Recall that the duration of the wave is qualitatively proportional to 1/ √ a; he duration of the 2nd wave in India was smaller than that of the 2nd-3rd waves in the USA and Europe.

The coefficient c coincides with that for the 3rd wave in the USA; it is much smaller than the one for the 1st wave in India, which was 5.75. Generally, this means that the population of India reacted faster (reduced the number of contacts) when the 2nd wave arrived. The selfimposed restrictions are very efficient hard measures. Biologically, the strain became more virulent. For instance, the likelihood increased dramatically for the whole family to become infected if one of its members is infectious vs. the 1st wave according to Indian medical officials.

We conclude that much greater number of detected cases is obviously related to the jump of the magnitude C. This is generally simply a scaling coefficient. However, when comparing different waves in the same country, it contains valuable information. Our theory results in the following output: between the waves, the virus in India acquired significantly greater ability to penetrate practically all strata and age groups of this very large and diverse country.

In the USA and Europe, the changes between the first 2 waves were mostly due to the relaxation of the hard measure (including selfrestrictions), i.e due to a very significant drop of the coefficient a. Some . CC-BY-NC-ND 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. increase of c was of the same origin. However in India, c diminished and a increases, so the virus itself evolved: became a "broader one". Generally, the changes of c, a, C can be used to analyze the trends in the virus's evolution. Though they are very much linked to our response to the threat in our modeling.

USA: waves 2 and 3. We use them heavily to adjust the forecasting component of our theory, especially the 3rd wave. The latter was used for (relatively) long-term forecasting, which appeared quite possible.

The 2nd wave in the USA: 06/16 -9/12, 2020. The two-phase solution worked well for the second wave in the USA. The accuracy is comparable with what we had above for the first waves in Japan, Italy, Germany, the Netherlands and UK. Upon subtracting 2.1M , the second phase matched well the following functions: (1, t) )), a = 0.06, d = 0.435.

We note that the initial transmission rate was c = 2.2 for the USA during the first wave. The parameters c, C and 0.6 in the first formula were determined for the period marked by red dots; the black dots form a control period. See Figure 11 . The projected saturation for u B . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 :

: for in the USA Figure 11 . 2-phase solution for the 2nd wave in the USA. is given by the formula t end = exp

) . Numerically, t end = 17.8463, which is 178 days from 06/16: December 11, 2020. Though this did not materialize since the USA entered the 3rd wave in the middle of September.

USA: the 3rd wave, 2020-21. This wave deserves a detailed analysis. It was used for long-term forecasting, which we will describe in consecutive steps, following the corresponding projections.

This wave was on top of the unfinished 2nd wave, so we subtract the starting total number of infections, which was about 6.9M on 9/24, when our curve begins. The red dots used to determine the parameters of u(t) were taken from 9/24 to 11/17/2020, when the 3rd wave in the USA still did not reach the turning point.

The 1st control period (black dots) was till 12/13. The match was already very good, See Fig. 12 . However the accuracy became even better in March. Possibly some early signs of herd immunity contributed.

The formula for u(t) for the 3rd wave was as follows:

The new c increased from c = 2.65 for the 2nd wave to 2.85 in a way similar to the passage from the 1st wave to the 2nd. The parameter a . CC-BY-NC-ND 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. significantly dropped from a = 0.06 for the 2nd wave to 0.02, 3-fold. We note that a = 0.06 was about 1/3rd of a = 0.2 for the 1st wave, so the same tendency persisted. Recall that the duration of phase (A) is qualitatively ∼ 1 √ a .

. CC-BY-NC-ND 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. Forecasting the saturation. According to our analysis in various countries the passage to the end of the current cycle of Covid-19 began in the middle of March, 2021. This is a stage when the levels of immunity and vaccination against the current strains become high enough. So it lasted for about 1.5 years with 2-4 waves in many countries. This is not in all countries, and the next cycle can be close-by.

The projection was 03/05/2021 for the saturation of the 3rd wave in the USA. It was supposed to be followed by some period of modest (essentially linear) growth under the (B)-mode (the second phase). This projection was posted on 11/17/2020; see Fig. 13 . The match during the first control period (black dots) was good.

The next control period is presented in Fig. 14. It shows that the accuracy of this forecast appeared even better in March. This is remarkable, because the previous waves in the USA were with some "unusually long" linear-type periods in the middle. Forecasting was more difficult with the USA than in almost any other countries we considered. Note that the massive vaccination began in the USA somewhat later, after the saturation on 03/05/2021: mostly from the end of March.

The 3rd wave in the USA was followed by some break until the 4th wave began. The 4th one was exceptionally short (for the USA). The country already reached the first stages of the herd immunity (for the current strains) in the middle of March, and the vaccination program . CC-BY-NC-ND 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. contributed very much, as it happened in Israel, UK, etc. The current cycle of this pandemic in the USA and in quite a few other countries is hopefully on its way to the end, though this is of course only for the strains of Covid-19 that dominate now and with various reservations.

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 Early projections. Let us provide some conditional projections for the latest waves in South Africa (ZAF ) and UK: Fig. 15, 16 . Both are obtained before the turning points, so the main parameters are not reliable. We take their "most plausible combinations"; generally, a "cone" of possibilities can be obtained.

The reasonable number of days necessary to make projections is of obvious importance. There were 37 days for the South Africa and just 22 for UK. The dots present the corresponding total numbers of detected cases diminished by that at day 1. The "forecasts" are in terms of the absolute total numbers of detected cases. The purpose is to test the limits of the usage of our theory as a forecasting tool.

It makes sense to provide the projection for the daily number of detected cases based on the usage of the R-number. Fig. 17 is from https://www.bbc.com/news/uk 57447632 (06/12/2021), with a reference to the UK Department of Health and Social Care. It gives the number of new daily detected infection greater than 15K at 06/21/2021. Our estimate is about 9.7K for the maximum of the new daily detections until the beginning of August, with the turning point around the beginning of July. According to Fig. 17 , the number of new daily infections will reach 40K in the beginning of July. Let us quote: "At the current rate of growth, as estimated on Thursday (06/10), the UK will reach 15,000 cases a day by 21 June and January levels of infections by late July -that is without any further relaxation of rules." Auto-forecasting. We mostly did this for the USA and Western Europe, but any countries can be "processed" during their second phases (for any waves); currently, there is no software for the 1st phases, i.e. that for the Bessel-type modeling.

The USA: from 03 to 05, 2020. We will provide the automated forecast for 50 states based on the period 03/17-05/27; the data were from https://github.com/nytimes/covid-19-data . Every state was processed individually with the interaction; see [Ch1] . Our approach to incorporating the interaction is of independent interest: we allow the curves for individual states to become decreasing as far as the total sum increases, which is motivated by physics.

Our program focuses on the last 20 days; however, the match with the total number of detected infections appeared perfect almost from 03/17/2020. See Figure 18 and compare it with Fig. 20 .

. CC-BY-NC-ND 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 sum of the a ically found curves , where the interaction s allowed. The red dots: the infections .

The B-type curves were found on the basis of the ed well . This means that the forecasting remained very stable for 2 months. Figure 18 . USA, the sum of the curves for individual states.

Such a high level of stability is actually rare in any forecasting, which made the expected chances good for reaching the saturation around 9/19/2020, which was the projection based on Fig. 18 . The saturation is of technical nature: it does not mean the end of the wave. Normally, it is followed by a period of modest linear growth of the total number of infections, which we model using u B (t). There are always remaining . CC-BY-NC-ND 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 June 22, 2021. Figure 19 . An auto-forecast for Europe as of 7/14 (2020) and new clusters of infection, though sometimes mode (A) alone is almost sufficient until the end of the wave. However, the hard measures were significantly reduced in the USA at the end of May practically in all 50 states. As a result, the number of states that reached phase 2 dropped from about 22 at 5/27 to 8 at 7/12 (2020). Then, in the second half of June, the USA entered the second wave. The program was quite stable for the 2nd wave in the USA ... before it entered the 3rd wave in the middle of September. We note that it was written for the 2nd phase only, when constant (daily) monitoring and tuning the process of management is more relevant. This is the period when the hard measures are being reduced.

Europe: summer 2020. The situation was quite stable in Europe in summer of 2020. We provide in Figure 19 a sample forecast our automated system produced for Western Europe till the end of July, 2020, to be exact, for the following 45 countries: Albania, Andorra, Austria, Belgium, Bosnia and Herzegovina, Bulgaria, Croatia, Cyprus, Czech Republic, Denmark, Estonia, Faeroe Islands, Finland, France, Germany, Gibraltar, Greece, Guernsey, Hungary, Iceland, Ireland, Isle of Man, Italy, Jersey, . CC-BY-NC-ND 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 June 22, 2021. Here and below the curve average is the maximum and the corresponding value of the average of 9 last curves u B (t) for "Europe". I.e. we consider the average of 9 curves we obtain for 9 consecutive days and then determine its maximum. The 9-day average is the simple average of the corresponding maxima for these 9 curves; i.e. it is the moving average. The main source of Covid-19 data we used was: https://ourworldindata.org/coronavirus.

As of July 8, 2020 the following countries had clear second phases: Albania, Bosnia and Herzegovina, Bulgaria, Croatia, Czech Republic, Greece, Kosovo, Luxembourg, Macedonia, Montenegro, Romania, Serbia, Slovakia, Slovenia. The forecasts were sufficiently stable, though Sweden, Poland, Portugal and some other countries did not reach phase 2 at that time. Such stability changed in fall 2020 due to the end of the vacation periods and the beginning of the school year.

Summer 2021. Fig. 20 provides an auto-generated projection for the USA as of 06/09/2021. This is without the state-by-state analysis: the total numbers of detected infections for the whole country are used. The projections are supposed to be constantly renewed during phase 2. The fact that the curve of (actual) detected cases was very well approximated by our u B (t) for the whole period of 90 days is very remarkable. Only the last 20 days are used for finding the (current) parameters. This is a strong confirmation of our two-phase solution. The USA was in the middle of the intensive vaccination program during this period. Also, it was already close to the initial levels of herd immunity (for the current strains) in March-May. The aggressive vaccination program and continued restrictions of contacts (though significantly relaxed in May-June) influenced this period.

For Europe (summer 2021), a similar projection is provided in Fig.  21 . Again, a perfect match with the actual data is for the whole period of 90 days; only the last 20 days are used to find the parameters.

With both forecasts, the growth of the number of the cases is expected this Fall due to the impact of the Delta-strain and other factors. This strain already influenced UK; see Fig. 16 . Europe and the USA are likely to follow. The Delta-strain might be a transitional one between cycles. Though some recombination of the existing strains . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 is more likely to become dominant. Really new strains can of course emerge, including those "responding" to the vaccination programs.

The situation with Covid-19 remains far from "normal" in many countries. However the levels of herd immunity reach sufficiently high levels and there are intensive vaccination programs in quite a few countries. These factors can be expected to end the current cycle, but this is of course conditional. Much will depend on the evolution of the strains of Covid-19 and there are already some new waves with a potential to trigger a new global cycle, which can be close-by.

However, a combination of hard measures, herd immunity, and the vaccination obviously began to work, though it took and will take enormous amount of efforts and a lot of time. Societal costs of the pandemic were already huge, and this is not the end.

Hopefully we will be better prepared to new cycles of Covid-19, including better mathematical tool for forecasting the spread of epidemics, understanding the uniformity of the waves of the infections, finding the key and meaningful parameters that determine the spread, and making solid mathematical predictions.

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 Figure 21 . A sample forecast in Europe: 06/08/2021.

Modeling Covid-19 appeared quite a challenge for existing mathematical methods, which are mostly based on the SID-type approach, suggested in the early 20th century. The following features of Covid-19 obviously require new methods.

(1) The curves of total numbers of detected infections are mostly of power-type for Covid-19, where the exponent diminishes over time.

(2) The saturation of the initial waves of Covid-19 is mostly because of the protective measures, not due to the herd immunity.

(3) The range and intensity of the protective measures used to fight Covid-19 are exceptional in the history of epidemics.

These factors are not really new in epidemics, but they do require a new mathematical theory. The prior approaches appeared insufficient for modeling the spread of Covid-19. For instance, the power-type growth of the total number of infections cannot be addressed within the SID-type models.

Our theory seems the first one when the power growth of the spread and the active epidemic management, including self-imposed restrictions, are considered the major factors. These assumptions result in differential equations depending only on the initial transmission rate . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 and the intensity of the protective measures. These parameters make perfect sense theoretically and practically, and can be measured reliably during relatively early stages of the waves of Covid-19.

The actual graphs of the total numbers of detected cases in many countries (all we considered) are described uniformly and with surprisingly high accuracy by our curves. The 2-phase solution is a combination of the Bessel-type curve for phase 1, which is the key in our approach, and its certain version for phase 2 (in terms of elementary functions). We note that the corresponding saturations were mostly due to the active management till summer 2021. They are of unstable nature: reducing the protective measures may result in the recurrence. Modeling and forecasting this kind of saturations requires sharp mathematical tools, which we hopefully provide.

Methods. The starting point of our approach to modeling the total number of infections during epidemics is the power growth hypothesis, which has solid confirmations with Covid-19, practically in all countries (apart of small initial periods). This is without considering epidemic management; we deduce it from the principle of local herd immunity.

The saturation of the corresponding waves of the spread of Covid-19 within a cycle is mainly due to the protective measures, the process of active management of the epidemic including our own actions. Protective measures are of course not unique for Covid-19, but their range and intensity reached unprecedented levels. Our model connects this kind of saturation with the asymptotic periodicity of Bessel functions, one of the deepest results in their theory. This is very different from the classical approaches based of SID, SIR, SIER models and some their variants, and the approaches in neighboring segments of ecology: invasion and interaction between species.

Due to a very limited number of parameters, actually 3 for our twophase solution, our model is much more rigid than any other ones. These parameters are quite meaningful mathematically and epidemiologically.

The key are c, the initial transmission rate, and a, the intensity of hard protective measures, which can be determined reliably at relatively early stages of particular waves of the epidemic.

These parameters are generally different for different waves, but we see some patterns here. The second waves are almost always for somewhat higher c and significantly lower a, the intensity of hard protective measures. The 2nd wave in India was exceptional: c dropped and a . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 increased. The number of contacts was reduced faster during the 2nd wave, much broader than the 1st one, and the protective measures, including important self-imposed restrictions, were significantly stronger. As for the 1st wave, the match with our u(t) was almost perfect.

We obtain a very good match practically for the whole periods of the waves of Covid-19 in many countries. It is actually surprising for such stochastic processes as epidemics, and of course this is a strong confirmation of our assumptions.

Since our theory was created and posted in the middle of April 2020, which was mostly in the middle of the 1st waves of Covid-19, we had a unique opportunity to test it in the course of epidemics. Namely, we systematically determined the parameters of our curves during relatively early stages and then tested the corresponding projections extensively for sufficiently long control periods.

Our usage of control periods is similar to routine testing the quality of the models used for forecasting share-prices in stock markets, where no approach can be accepted without real-time runs and carefully crafted historic experiments that exclude any "usage of future" as far as possible. This kind of "discipline" is not present in forecasting the epidemics, at least by now. The results of checking our theory during the control periods, including automated forecasting programs we developed, are an important part of papers [Ch1, Ch3] . An example of the long term forecast is the 3rd wave in the USA discussed above.

We demonstrated that Bessel-type functions describe very well the growth of the total number of detected cases in many countries. Mathematically, we successfully model the passage from ∼ t c , describing the initial growth of the total number of detected cases, to ∼ t near the turning point, and then almost all the way to the saturation of the current wave, when c drops almost to 0. Here c is the initial transmission rate, which can be captured at relatively early stages of the waves of Covid-19.

Our functions are solutions of ODE describing the epidemics under active management, the system of protective measures. The selfimposed restrictions on the contacts during the epidemic are and always were important protective measures. Constant receiving reliable information on the course of the epidemic by the population is of obvious significance here. The travel restrictions, closing the places where . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 the spread is the most likely, and the vaccination are classical state controlled measures.

The saturation due to active protection measures is of unstable nature, unless the herd immunity becomes a significant factor. Its forecasting requires an exact mathematical theory, which we try to provide. There will be an endless discussion of the efficiency of different measures and different management approaches until verifiable trustworthy mathematical models and the corresponding software are developed and implemented practically.

The verification of any models, including this one, does require algorithms that can be used by anyone, not only by their creators, the ultimate test of their validity. This is one of the reasons why we wrote our programs. They are posted in [Ch1, Ch3] and can be used by anyone to model any countries and regions, though only for the late stages of Covid-19 so far (mode (B)).

The new theory seems a solid basis for reaching the next level, which is forecasting. It already describes the curves of total numbers of detected infections with high accuracy and with surprisingly high level of forecasting and stability of the auto-projections for later phases. Though forecasting is always a challenge.

Basic parameters. The small number of the parameter we employ explains well the uniformity of the curves of total numbers of detected infections of Covid-19 in many so different countries, as well as the mathematical similarity of different waves. They are:

(1) the initial transmission rate c, which can be determined at relatively early stages of the current wave, (2) the intensity of hard measures a, which become sufficiently stable near and after the turning point, (3) the intensity b of the measures during the second phase, which begins near the end of phase 1.

In contract to c, the intensity of the measures is of course more time-dependent, but a appeared sufficiently stable for long periods. Concerning b, it must be adjusted constantly at the later stages, because no country is really isolated and the management of the epidemic, including the self-imposed protection measures, becomes more flexible at these stages, so less predictable mathematically.

The scaling coefficient C of u(t) is adjusted to match the real numbers of cases. The coefficient of u 2 , the non-dominant solution of our . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 system of ODE, is another parameter. It is mainly responsible for the "effects of the second order" and does not seem really critical for forecasting; the dominant Bessel-type solution u 1 is expected to be sufficient for this. There are technical reasons why C increases when u 2 is added with a significant positive coefficient: u 2 (t) becomes negative for sufficiently large t and contributes greatly for such t.

The coefficient C is generally a technicality, to make the output convenient to present as a graph. However it provides valuable information when different waves of the epidemic are compared in the same country. An example is our analysis of the 1st and 2nd waves in India. The whole country was affected during the 2nd wave, when only some strata of this huge and diverse country were mostly infected during the 1st wave. The coefficient C dramatically increased during the 2nd wave there, when the coefficients c, a became similar to those in Europe during the 1st-2nd waves. The match of our u(t) and the total number of detected cases was very good with the 2nd waves.

Mathematically, u 2 to captures some features of the curves of total numbers of detected infections, especially during the initial stages. It can be "seen" in quite a few countries, though u 2 mostly influences the effects of second order. From the view point of our approach, u 1 and u 2 , both, could be expected to occur, which was really confirmed. See [Ch1] for some discussion.

The fact that we were able to describe such complex stochastic processes as epidemics with only 3 basic parameters seems a real discovery. This worked very well almost everywhere (in all countries we considered), but it will take time to understand the scope of this new theory and to begin using it practically. Obviously such a remarkable universality can have deep roots in general biology and the theory of random processes. The same system of ODE serves quite a different situations, including momentum trading in stock markets. It models MRT, momentum risk-taking. Its "deduction" is different for epidemics (and in invasion ecology). The price function in momentum trading is replaced by the protection function p(t) for (1).

Beyond this paper. Needless to say that the vaccinations makes us closer to the herd immunity, at least for the current strains of Covid-19. There are many challenges here [AVTC] , and, anyway, the control of the efficiency of the vaccination programs requires exact mathematical . CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 tools; there are other factors. In our approach, this means measuring their impact on the parameters, especially on c, the transmission rate.

The current waves of Covid-19 can be hopefully the last for this particular cycle, assuming some stabilization of the strains involved. The forecast we posted in November 2020 for the USA and presented in Fig. 13 was with the expected saturation (of phase 1) on March 5, 2021. Its accuracy appeared very good. The total number of all cases, detected or not, is of course greater; the impact of herd immunity was and can be expected at such levels of detected cases. The vaccination did not play a significant role in the beginning of March, but certainly it does now (Summer, 2021) . This is a major protective measure; however its success requires forecasting future strains, which is difficult.

Concerning the new strains, presumably about 6 months are necessary for them to reach a significant presence counting from their "first appearance". The latter can be the moment of their creation, which can be expected at late stages of the prior waves, or when new strains are "imported" from other regions. Also, the immunity upon the recovery from Covid-19 or upon the vaccination is expected to begin to fade after about 6 months, which is preliminary too.

Covid-19 demonstrated a robust potential for creating various VOC, variants of concern. Its ability to evolve resisting the vaccines and the treatment remains to be seen, but it can be expected to be ample. Generally, it can take only several infected individuals where the virus can "stay" for 1-3 months to generated many mutations, but fortunately very few of them are stable and have a chance to reach any dominance. The protective measures greatly restrict their spread (always did); a mathematical theory of their impact combined with power-law of epidemics is the purpose of [Ch1, Ch3] and this paper.

Basically, the end of a cycle is when the current strains become less infectious due to growing immunity, vaccinations, the change of the season and similar factors. The epidemic can continue between the cycles in other places, zones, and zoonotic reservoirs, which generally results in a significant renewal of the strains. We hopefully approach the end of the 1st cycle in the Northern Hemisphere within the current dominant strains due to a combination of herd immunity, vaccinations and summer time, though this is far from certain. The biological aspects are of obvious importance here, which are beyond this paper.

. CC-BY-NC-ND 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 June 22, 2021. ; https://doi.org/10.1101 https://doi.org/10. /2021 Mathematically, a challenge is to understand the formation of waves within one cycle. Recall that we stop the usage of Bessel functions when the corresponding u(t) begins to diminish, which is impossible for the total number of cases. Then we switch to mode (B). Given a country, the waves follow some patterns. However the pauses between them and the parameters c, a, C fluctuate significantly. It seems more predictable when the next wave is due to the relaxation of the protective measures at the end of the previous one, but there can be other factors, including new strains. Anyway, 2-4 waves (in one country) are insufficient for modeling their recurrence, theoretically and practically.

Main practical applications. They seem for modeling individual waves. Namely, the natural aims are as follow:

Creating software for monitoring the dynamic of the epidemic. This is especially needed at the later stages of its waves (2nd phases in our terminology), when tuning the management becomes important. (3) When comparing different waves in the same country, conclusions can be made about the changes in the epidemic management (including self-imposed restrictions), and the changes of the virus behavior. (4) Extending this theory from modeling epidemics to Invasion Ecology, more specifically, to transient processes of the invasion and interaction between 2-3 species; applying it to momentum risk-taking. Acknowledgements. The author thanks very much David Kazhdan for valuable comments and suggestions; a good portion of this paper presents author's attempts to answer his questions. Many thanks to Eric Opdam and Alexei Borodin for their kind interest and ETH-ITS for outstanding hospitality

Author's Affiliation

Challenges in creating herd immunity to SARS-CoV-2 infection by mass vaccination

Least-squares finite element method for a meso-scale model of the spread of Covid-19

A mathematical model reveals the influence of population heterogeneity on herd immunity to SARS-CoV-2

Are RNA viruses candidate agents for the next global pandemic? A review

First passage events in biological systems with non-exponential inter-event times

Momentum managing epidemic spread and Bessel functions

Artificial intelligence approach to momentum risk-taking

A surprising formula for the spread of Covid-19 under aggressive management

Modeling infectious disease dynamics

The many guises of R 0 (a didactic note)

Mathematical Tools for Understanding Infectious Disease Dynamics, Princeton Series in Theoretical and Computational Biology

Factors that make an infectious disease outbreak controllable

Transients: the key to long-term ecological understanding?

The mathematics of infectious diseases

Periodicity in Epidemiological Models

Bayesian-based predictions of COVID-19 evolution in Texas using multispecies mixture-theoretic continuum models

Thinking, fast and slow

Memorybased meso-scale modeling of Covid-19-County-resolved timelines in Germany

The Mathematics Behind Biological Invasions

Strong correlations between power-law growth of COVID-19 in four continents and the inefficiency of soft quarantine strategies

Power-Law models for infectious disease spread

Graph theory suggests COVID-19 might be a "small world" after all

Epidemic psychology: a model

A network-based explanation of why most COVID-19 infection curves are linear

Simulating the spread of COVID-19 via a spatiallyresolved susceptible-exposed-infected-recovered-deceased (SEIRD) model with heterogeneous diffusion

A Treatise on the Theory of Bessel Functions