key: cord-0654684-x1jtn0px
authors: Lukyanets, S. P.; Gandzha, I. S.; Kliushnychenko, O. V.
title: Modeling and Controlling the Spread of Epidemic with Various Social and Economic Scenarios
date: 2020-06-12
journal: nan
DOI: nan
sha: 21e2d8e35d6a32853f2287ae28d8b01f3eade2af
doc_id: 654684
cord_uid: x1jtn0px

We propose a novel model for describing the spreading processes, in particular, epidemics. Our model is an extension of the SIQR (susceptible-infected-quarantine-recovered) and SIRP (susceptible-infected-recovered-pathogen) models used earlier to describe various scenarios of epidemic spread. As compared to the basic SIR model, our model takes into account two possible routes of virus transmission: direct from the infected compartment to the susceptible compartment and indirect via some intermediate medium or fomites. The transmission rates are estimated in terms of the average distances between the individuals in selected social environments and characteristic relaxation times. We also introduce a resource activation function that reflects the load of the epidemics on economics and the limited capacity of the medical infrastructure. Our model brings an advantage of building various control strategies to minimize the effect of the epidemic and can be applied to modeling the recent COVID-19 outbreak.

As widely acknowledged [1] [2] [3] [4] [5] , the spread of information, rumors, ideas, or deceases share many similarities and, in most cases, spreading processes are governed by similar models.

There are entire classes of models in terms of complexity or approximations to describe the spreading processes, starting from dynamic systems, e.g. [6] [7] [8] [9] [10] [11] [12] [13] , that stem from classical works [14] [15] [16] [17] [18] . More complex models include stochastic effects, e.g. [19] ; models with spatial flows, e.g. diffusion [20] ; and models allowing for nontrivial spatial structure or topology, e.g. [1, 21] . One common feature in the majority of these models is the presence of kinetic coefficients or parameters that characterize the probability of elementary processes or reactions per unit time. On the one hand, these model parameters (constants) determine the instability points and characteristic rates of instability growth. On the other hand, however, the kinetic coefficients or probabilities are not the proper control parameters that are readily available, and they can only be estimated indirectly in terms of other parameters.

The dynamics of any spreading process, e.g. an epidemic, is determined by the individual peculiarities of people in the selected social group, by the type and mechanism of infection, etc. To describe such a process, it is necessary to indicate some ways for a specific influence on the process by selecting, in particular, a proper quarantine strategy. A dramatic example is the spread of COVID-19, when different countries resort to different quarantine measures [22] . Such measures are controlled by the choice of certain parameters available to society (e.g. working day duration, average density of people in public places, frequency of disinfection, intensity of * lukyan@iop.kiev.uagandzha@iop.kiev.uakliushnychenko@iop.kiev.ua transport communications, etc.). The problem of strategy selection reduces to problems of optimal control theory for feedback systems (theory of games in the more general case [23] ). An important point in strategy selection is to single out the control parameters that are readily available in the system.

One of the basic constants in spreading models is the reaction cross-section or its characteristic rate. For instance, the infection reaction is associated with scattering of the infected individuals on the spatiotemporal fluctuations of population density. As is well-known [24] , the infection probability depends on population density non-monotonically. At low densities, the low probability of scattering of the infected individuals on the susceptible individuals results in a small total cross-section of the infection reaction. At high densities, the mobility of spreaders decreases, and it mitigates the spreading process. Moreover, the scattering of spreaders on local population-density fluctuations should determine the critical concentration of secondary infected cases sufficient for the initiation of the collective process, i.e. epidemic. In general, kinetic description of spreading is quite a complex problem that needs to account for the changes of the internal state of the spreaders in the course of collisions as well as the presence of spatial inhomogeneity in the system and its non-equilibrium properties.

parameter E serving as its height. Therefore, it can be assumed that the rehabilitation rate should have an activation character, ∼ exp(−E/ρ), similar to the temperature dependence of the activation process with activation energy E, see, e.g., [26] [27] [28] .

The presence of activation process indicates that the system can be subjected to the so-called explosive instability. For example, when a chemical reaction occurs with the release of heat and has an activation character, it goes faster at higher temperatures. This leads to yet greater temperature increase and, as a result, to the so-called Zeldovich-Frank-Kamenetskii explosive instability [29] [30] [31] . The fight against the epidemic involves additional costs, e.g. associated with quarantine measures, which result in the reduction of the production of the collective formal resource ρ. The resource being depleted, the quality of medical services and the rehabilitation (recovery) rate drop. As a result, the number of active members in the population goes down. This, in turn, leads to a further decline in the collective resource production, with the level of income needed for basic survival being lower and lower. Such a scenario finally results in the total collapse if the system.

In this paper, we try to address the points mentioned above. In other words, we are interested not only in modeling the spreading process but also in identifying the ways it can be controlled by socially available parameters and which consequences such control measures can lead to. To demonstrate a number of possible effects, we turn to the simplest feedback model that is based on dynamical systems. We would like to mention two features of this model.

(I) In order to take into account the fact that an individual can stay in different environments or public places (e.g. transport, shops, work) that are characterized by different densities of surrounding individuals, it is convenient to consider the selected social group and its average daily cycle (T ). For instance, we can introduce the average daily time slots the individual spends at home (T 1 ), at shop or other public places (T 2 ), in transport (T 3 ), at work (T 4 ), etc. and define (at least under normal conditions) the characteristic density c j of individuals or the average distance ℓ j between them. We also make a qualitative assessment of the dependence of the reaction (infection) rate on the density of individuals or average distance between them in order to clearly identify the set of control parameters {T j , ℓ j }. Finally, we take into consideration a possibility of indirect transmission through the medium (e.g. contaminated water [32, 33] ) or via the so-called fomites [34, 35] . The indirect transmission rate is determined by c j (ℓ j ) as well.

(II) The relationship between economics and epidemic is taken into account by means of the above-mentioned activation mechanism ∼ exp(−E/ρ) for the characteristic relaxation times. To this end, we use the simplest form of the equation for the dynamics of the formal resource ρ.

Before we proceed to the consideration of our main model described in Sect. IV, we make an attempt to understand the effect of the relaxation activation mechanism when the system can pass through a number of quasi-stable states and then collapse. In this case, its dynamics resembles the so-called devil's staircase [36] providing a clear illustration to the global pandemic scenario of the world never being the same again after the epidemic. To this end, we turn to a simple toy model that accounts for a formal demographic and economic resource and is more suitable for describing a primitive community or family.

We first consider a simplified model where some social group of individuals is subdivided into two compartments: susceptible (S) and infected (I). The susceptible individuals are infected at some transmission rate β, which is defined as a product of the contact rate and the probability that a contact of infected individual with a susceptible individual results in transmission. The infected individuals recover and become susceptible again with some recovery rate Γ i , which identifies the probability of recovery per unit time and is estimated as the reciprocal of the mean time spent in the infectious class. The recovery process is governed by the general economic situation described by some integral activation parameter E which reflects the cost of medical and other essential life services as well as the bare subsistence level of consumption. The corresponding mathematical model is given by the following two ODEs:

where the operator ∂ t stands for the derivative with respect to time t. Here s is the number density of susceptible individuals and i = 1 − s is the number density of infected individuals. The total number of susceptible and infected individuals is assumed to be constant. The initial conditions at t = 0 are taken as

i 0 being the initial number density of infected individuals. The function ρ represents some general "resource". Formally it corresponds to the collective product or earned money. The production of this resource per unit time is proportional to the number density of working individuals, s (the infected individuals are assumed to be not working). The constant G formalizes the resource volume generated by them per unit time. The second term, Γ ρ ρ, formally describes the collective expenses or taxes. Roughly speaking, the expenses are assumed to be proportional to earnings. Thus, the coefficient Γ ρ represents the resource consumption rate. The last term Λ represents the fixed expenses necessary for keeping some infrastructure (e.g. amortization, municipal services, etc.). In the case of unlimited resource (E ≪ ρ), the equation for s reduces to the basic SIS (susceptible-infectedsusceptible) model, whose solutions are well studied [10, 11, 17] . The purpose of this paper is to investigate the effect of nonzero activation parameter E (which we will also refer to as "activation energy") on the dynamics of epidemic described by system (1).

As we understand, if there are no black swan events like the epidemic, there exists a quasi-stationary equilibrium state with ρ = ρ (0) = const. It is called the disease-free equilibrium and is given by the trivial stationary solution of Eqs. (1):

On the other hand, under the stress situations like epidemic, the system can go out from the disease-free equilibrium, with resource depleted. Indeed, another stationary solution to the equation for ρ is given by

where s * is given by the following transcendental equation:

The parameter

is the basic reproduction number. It defines the average number of transmissions one infected individual makes in the entire susceptible compartment during the entire time of being infected. When R 0 1, the disease-free equilibrium is stable, and there is no epidemic outbreak. When R 0 > 1, the disease-free equilibrium is unstable, and the system evolves to the new equilibrium state {s * , ρ * } called the endemic equilibrium [10, 11] . When E = 0, Eq. (5) has two solutions. The first solution,

is stable at R 0 > 1 and defines the endemic equilibrium point. The second solution s * = Λ G is unstable for all R 0 . To simplify our further analysis for the case E > 0, we put Λ = 0. Then Eq. (5) can be rewritten in a simpler form:

where

For 0 < E < E c , where E c = e −1 ≈ 0.368, Eq. (8) has two solutions: z 1 > E c (which defines the endemic equilibrium point) and 0 < z 2 < E c (which is always unstable). For E = E c , there is one solution z 1,2 = E c . Finally, there are no real solutions for E > E c . The last case is most important for our consideration. It means that at some E and R 0 there is no stable endemic equilibrium because of resource depletion, and the dynamical system described by Eqs. (1) should collapse (such a situation is sometimes called the explosive instability). The critical activation energy at which the system starts to collapse is given by a formula

Except for the condition 0 E < E c , the endemic equilibrium point should also meet the requirement of s * < 1, which effectively implies that E < E e , where

Relations (10) and (11) define two critical curves in the (R 0 , E) plane which determine the evolution scenario for dynamical system (1) . Depending on the values of parameters R 0 and E, the system can evolve into three possible states: disease-free equilibrium, endemic equilibrium, or collapse. Figure 1 shows the corresponding phase diagram.

The above analysis is supported by the results of numerical integration of Eqs. (1) demonstrated in Fig. 2 (the case R 0 < 1) and Fig. 3 (the case R 0 > 1). In these examples we normalized the resource function by its trivial stationary value ρ (0) , which is equivalent to taking ρ (0) = 1. In this case we have ρ * = s * .

When E = 0, the dynamics of system (1) follows the basic SIS model. It evolves to the state of disease-free equilibrium at R 0 1 [ Fig. 2 (a)] and to the state of endemic equilibrium at R 0 > 1 [ Fig. 3(a) ].

When E > 0, some part of the resource is consumed, and the number of healthy (susceptible) individuals decreases [ Fig. 3(b) ]. There is a critical value E e defined by formula (11) at which the system evolves to the endemic equilibrium even at R 0 < 1 [ Fig. 2 (b) ]. This scenario is impossible in the basic SIS model. When the activation energy is further increased and becomes larger than the critical value E c defined by formula (10), the endemic equilibrium is no longer stable and the system collapses to the state s * = ρ * = 0 [ Fig. 2(c,d) ]. This means that all the individuals become infected and there is no resource to reverse the epidemic back. The same scenario is observed in the case

When E is above the critical value E c but still close to it, the system first tries to occupy the quasi-stationary endemic state (which is unstable). This process can take quite a long time and then the system finally collapses The simplified SIS-like model and example considered in this Section demonstrate that in the case of limited resource (E > 0), there exists a certain critical point for any basic reproduction number R 0 at which the system collapses and can no longer stabilize and return to the stable pre-epidemic or after-epidemic state. This fact provides a clear illustration to the global pandemic scenario of the world never being the same again after the epidemic.

As mentioned above, the spreading process can be controlled only by accessible for society parameters. These parameters are, e.g. the mean distance between individuals and the mean time spent in a particular place.

In order to take into account the fact that an individual can stay in different environments or public places (e.g. transport, shops, work) that are characterized by different densities of surrounding individuals, it is convenient to consider the selected social group and its average daily cycle (T ). For instance, we can introduce the average daily time slots the individual spends at home (T 1 ), at shop or other public places (T 2 ), in transport (T 3 ), at work (T 4 ), etc. and define (at least under normal conditions) the characteristic density c j of individuals or the average distance ℓ j between them. We also make a qualitative assessment of the dependence of the reaction (infection) rate on the density of individuals or average distance between them in order to clearly identify the set of control parameters {T j , ℓ j }. We estimate the reaction cross-section, using a simple approach, in order to take into account a safe distance, and consider also the additional mechanical infection associated with infecting by the environment. The later is also depends on the population density in a given place. Such approach allows not only to single out the characteristic parameters one can influence, but also to describe, in the same respect, another social groups, e.g. the category of people that go to work in the subway, with high population density, or by car, taking into account their intersecting points. When determining the chemical reaction cross-section via a field or some mediator, we deal with parameters which should be obtained from general considerations. For example, the probability to be infected per unit time for an individual in the vicinity of infected one for a certain period of time which depends on the individual organism resistance.

Here, we roughly estimate the dependence of mean infection rate ω j (reaction rate) on characteristic distance between people for an individual in a "social" region (location) "j" with a given mean population density c j ∼ (ℓ j ) −2 , where ℓ j is the mean distance between them.

In general case, the probability of being infected ("interaction potential") depends on the distance between infected and susceptible, and, as may be assumed, is a damped function with a characteristic correlation radius (decay length) ℓ c . Then, the scattering (reaction) crosssection is determined by the probability of finding the infected individual at a certain distance from susceptible one and by the interaction potential. By averaging over all the configurations of spatial locations of infected, susceptible and immune individuals, we obtain the characteristic infection rate ω j which depends on ℓ j and ℓ c . In the case of non-trivial spatial structure (topology) of population, the scattering problem is correctly described by the two-point correlation function or so-called structure factor.

However, for simplicity, here we assume that the elementary infection act probability for susceptible individual, falling into the vicinity of radius ℓ c around the infected one, does not depend on the distance and equals ω. We also assume that the event of being infected from two different infected individuals are independent. In what follows, ℓ c is associated with socially safe distance.

For simplicity, we suppose that, for each social zone j, the mean population density is fixed under any social conditions. The mean number of persons falling into the area of radius ℓ c around susceptible individual is N j ≈ c j πℓ 2 c . One can estimate the infection probability of susceptible individual per unit time ∆P j = ωq(N j − 1), where q ≈ I/N is the probability to be infected. By summing over all the susceptible individuals in the system, we obtain the characteristic growth rate of the number of infected individuals, S∆P j , in the system with population density c j . By dividing on the total population number N , the rate of infection reaction in the region with density c j takes the form

i.e.

Here s is the number density of susceptible individuals and i is the number density of infected individuals. We can estimate the relative error for ω j caused by fluctuation of the number of individuals falling into the area of radius ℓ c by using, for simplicity, a triangular lattice with a lattice constant ℓ 0 (characteristic size of the individual). The deviation of squared number of particles

where θ is the probability that lattice site is occupied by a person, and N m is the maximal number of persons that can fall into the area of radius ℓ c :

Then, the relative error for ω j can be estimated as

assuming that ℓ c ∼ 4m, and minimal distance between persons ℓ 0 ∼ 1m, and ℓ 0 < ℓ j < ℓ c . Generally, it should noted that characteristic reaction cross-section can have more complex dependence on density ∼ c α , 1 ≤ α ≤ 3, and other type of non-linearity than si [57].

Next we roughly estimate the probability of elementary infection act ω. Let ℓ 0 be the minimum possible distance between two individuals. Assume that the suspected individual is infected with probability p 0 = 1 if he stays at the distance ℓ 0 from the infected individual for some characteristic time τ is (e.g. 8 hours). Then the probability of becoming infected at the minimum distance per unit time is ω = a/τ is . The parameter a can be estimated as the ratio of the "contact" area between the susceptible and infected individual to the total area determined by the correlation radius ℓ c : a ≈ ℓ 0 /ℓ 2 c . Then we assume that each individual can stay in N ω different locations during some typical period of time T (e.g. one day). Each location is characterized by its own average distance ℓ j . Then the integral transmission rate β for all locations is given by a formula:

where T j are the characteristic times spent in each of the locations j, with

In practice, we will limit ourselves by four locations (N ω = 4). Location 1 refers to staying at home (limited social contact), location 2 refers to shopping and other social contacts during a day, location 3 refers to staying in public transport (where the transmission probability is the highest), location 4 refers to staying at work.

Our estimates of the transmission frequencies ω j and rates β for various sets of social control parameters T j and ℓ j are given in Table II (see Appendix A) . We considered three possible scenarios: casual (no restrictions), soft quarantine (social distancing in place), and strict quarantine (restricted public transport, limited social contacts, partial cutting of economic activity). The soft quarantine measures reduce the integral transmission rate β by a factor of 3, and the strict quarantine measures reduce it by a factor of 7. Thus, our rough estimates given by formulas (12) and (14) allow the transmission rates to be controlled by the proper choice of parameters T j and ℓ j in different social scenarios. Our β estimates in the case of casual scenario fall in the range derived from the statistical data for the COVID-19 epidemic in the Wuhan city [37, 38] .

Here we consider a more realistic model of epidemic evolution. In addition to the direct transmission of infection by direct contact of the susceptible individual with the infected individual, we also consider the indirect transmission of the pathogen (e.g. virus) either through the medium (e.g. contaminated water [32, 33] ) or via intermediate objects (like hands, counters, doorknobs, etc.) often called fomites [34, 35] . Such an intermediate medium or fomites will further be referred to as "cloud". We associate a separate cloud with each social location j described earlier in Section II.

We subdivide the selected social group into five compartments: susceptible (S), infected (I), quarantine (Q), recovered (R), and deceased (F). The susceptible individuals (S) are those who are at risk to be infected. The infected individuals (I) are those who have been infected and pass through the virus incubation period. Then they either recover (with or without the acquired immunity) or pass to the severe form of decease (like fever or organ disfunction) when they need to stay isolated either at home or at hospital getting a medical help. Such individuals are referred to as the quarantine individuals (Q). The quarantine individuals either recover (with or without the acquired immunity) or die. The corresponding mathematical model is given by the following system of ODEs:

where the functions s(t), i(t), q(t), and r(t) describe the number densities of the susceptible, infected, quarantine, and recovered compartments, respectively. Here the index j = 1, . . . , N ω runs through all the social locations, and N ω is the total number of such locations taken into consideration (see Section II). The average frequencies ω j of infection transmission in each of the locations j are defined by formula (12) . The typical characteristic times T j spent in each of the locations j are given in Table II , the total time spent in all the locations being constant each day and equal to the day duration T [see Eq. (15) ]. The number density of fatal cases is described by the function f (t) that is determined by the equation

Each of the functions d j (t) describes the density of pathogen per individual in the cloud j and contributes to the indirect transmission of the infection from the cloud to susceptible individuals. The pathogen dynamics in the cloud j is given by the following equation:

(16c)

The first term, σ j i, describes the pathogen shedding by the infected individuals into the cloud, and the second term, γ j d j , describes the pathogen decay in the cloud due to natural inactivation, decontamination, or other routes. The resource function ρ and the activation parameter E have the same meaning as in Section III. The corresponding resource balance equation is

The above equations are supplemented with the following initial conditions

The total number density is assumed to be constant:

Effectively, our model is a combination or extension of simpler models like SIR (susceptible-infectedrecovered) [6, 7, 10, 11, 16] , SIRD (susceptibleinfected-recovered-deceased) [39] , SIQR (susceptibleinfected-quarantine-recovered) [40] , SIWR/SIVR/SIRP (susceptible-infected-recovered-pathogen) [32, [41] [42] [43] , and EITS (environmental infection transmission system) [34] .

B. Assumptions 1. The system is closed (no migration outside the selected population group and no one is added to the population).

2. The natural demography is ignored.

3. The uniform mixing of people is assumed (no uneven spatial distribution). The distribution of the pathogen in each cloud is assumed to be uniform. 4 . The latent period from exposure to the onset of infectiousness is ignored, i.e. all the exposed individuals are assumed to be infected and can infect others. The SEIR (susceptible-exposed-infected-recovered) model was demonstrated to have no practical advantage as compared to the SIR model [37] .

The full description and indicative values of the parameters of Eqs. (16) are given in Table III (Appendix A) .

In particular, the parameter Γ is describes a rate (we call it the relaxation frequency) at which the infected individuals recover without acquiring the immunity and come back to the susceptible compartment. The parameter Γ iq describes a rate at which the infected individuals develop the severe condition and pass to the quarantine compartment, where they become isolated either at home or at hospital. The parameter Γ ir describes a rate at which the infected individuals recover without complications and acquire the immunity. It is proportional to the probability µ of acquiring the immunity (see Table III ). The relaxation frequency Γ i is defined as follows:

where τ i is the characteristic pathogen incubation period. The parameters Γ qr and Γ qs describe the rates at which the quarantine individuals recover with or without the acquired immunity. The parameter Γ qf describes the fatality rate in the case of unlimited resource (E ≪ ρ). This case models a situation when the quarantine individuals all get the necessary medical care and the fatalities are only attributed to insuperable health complications (such as concomitant diseases or age factor). In the opposite case, when E ≫ ρ, the fatality rate is at its maximum and is equal to the total relaxation frequency for quarantine individuals,

where τ q is the average quarantine/hospitalization period. This case refers to the collapse of the medical system when the quarantine individuals get no medical help or treatment. Each of the parameters Ω j describes the typical frequency of the pathogen transmission from the cloud j to the susceptible individual. Two other cloud-related parameters, σ j and γ j , are the pathogen shedding rate (infected-to-cloud) and decay rate in the cloud, respectively. Their estimates are provided in Appendix A.

Finally, the constant G formalizes the resource volume generated by working individuals per unit time. It is proportional to the average working time T 4 for the given social group (see Table II ). The parameter Γ ρ represents the resource consumption rate. The parameter Λ represents the fixed expenses necessary for keeping some infrastructure (e.g. amortization, municipal services, etc.).

Note that the instantaneous number density of quarantine (ill) individuals q(t) is not often a convenient indicator for practical applications. The integrate number density of quarantine individuals for a certain period of time can be used instead. It is calculated as follows

In the case of unlimited resource (E = 0), the integrate number density of quarantine individuals in the end of epidemic (t → ∞) is proportional to the final number density of fatal cases,

where η is the probability of the fatal scenario for the quarantine individual.

For our further analysis, we first find the stationary solution to Eq. (16c) for the cloud j:

This relation allows us to introduce the dimensionless pathogen concentration in the cloud j, namely

so that p * j = i * . Then Eqs. (16a) and (16c) can be rewritten as

where β is the direct transmission rate (infected-tosusceptible) given by formula (14) and β j are the scaled indirect transmission rates via each of the clouds j (infected-cloud-susceptible):

The scaled indirect transmission frequencies ν j are defined as

Equations (22) supplemented with Eqs. (16b) and (16d) possess two equilibrium points. The first one is the disease-free equilibrium

with ρ (0) given by formula (3) . The second one is the endemic equilibrium

where the basic reproduction number R 0 is defined as 

is the aggregate indirect pathogen transmission rate via all the clouds. The parameter R 0 represents the number of secondary cases generated if a single infected individual is introduced into a completely susceptible population. When R 0 1, the disease-free equilibrium is stable, and there is no epidemic outbreak.

When R 0 > 1, the disease-free equilibrium is unstable, and the system evolves to the state of endemic equilibrium.

The stationary points for the resource function, number density of recovered individuals, and the number den-sity of fatal cases are given by the formulas

The variable resource is seen to have a profound effect on the endemic number densities of the recovered individuals and fatal cases. However, in contrast to the example considered in Section II based on a simpler SIS-like model, our extended model is more stable in regard to the resource depletion and does not allow for the total collapse scenario with s and ρ vanishing to zero at E ≫ ρ.

The first two equations of system (22) are strongly nonlinear. There are no analytical solutions known for the general form of these equations. However, in one particular case, when Γ is = Γ qs = 0 (µ = 1, no loss of immunity) and E = 0 (unlimited resource), one can get the following asymptotics at t → ∞:

Equation (30) can be used to control the accuracy of the numerical integration of system (22) . In the examples considered in the next Section, we used the fourth-order Runge-Kutta method to integrate Eqs. (22) numerically with step ∆t = T /10 sufficient to achieve the reasonable accuracy. In the case µ = 1 (no loss of immunity) our numerical estimate of s ∞ coincided with the value given by Eq. (30) to an accuracy of 10 −10 .

Now we consider particular examples to demonstrate various effects described by our model. First we focus on the case when there is no resource depletion (ρ ≫ E), so that the resource activation mechanism could be ignored (E = 0). We illustrate the "patient zero" phenomenon, then demonstrate the effect of the cloud, and finally model a quarantine scenario. Next we consider an example of a social group with limited resource (E = 0). Figures 4(a,b) show the number densities of the susceptible, infected, quarantine, and recovered individuals as functions of time in the case of two different initial number densities of infected individuals i 0 for the fixed basic reproduction number (R 0 = 2.5). Panel (a) corresponds to an initial density of one per 10 million, and panel (b) corresponds to an initial density of one per 100 thousand. The number density of infected individuals exhibits a typical peak and then drops. The peak has the same height for the both initial densities, but in the case of larger i 0 it is reached much faster (see Table I ). This example serves as an illustration of the "patient zero" phenomenon. When R 0 > 1, the epidemic spreads even when it starts only from one infected individual (patient zero). Then it reaches the same intensity in a certain period of time, which is shorter when the initial number of infected individuals is larger. This effect is also clearly seen in the phase portraits i(s) at different i 0 (Fig. 5) . Table IV (Appendix A) , and all other model parameters were selected according to Table III therein.

The corresponding parameters of the peaks are given in Table I . The inclusion of the cloud increased the basic reproduction number R 0 , so that the peaks of infected (i max ) and quarantine (q max ) densities become larger and shift to shorter times. The total number density of quarantine individuals (q Σ ) ∞ is significantly larger in the case with cloud.

The epidemic dynamics is governed by the basic reproduction number R 0 . The epidemic starts to spread when R 0 > 1. The greater R 0 , the higher are i max , q max , and (q Σ ) ∞ . Thus, to reduce the epidemic peak and to slow the epidemy down, one should reduce R 0 . According to formula (27) , this can be achieved by reducing the transmission frequencies ω j and, therefore, the transmission where t 1 is the moment when the quarantine starts, t 2 is the moment when the quarantine ends, β (0) and β (0) p are the direct and indirect transmission rates during the "no quarantine" period, and β (1) and β (1) p are the transmission rates during the quarantine period. Figure 6 presents the results of our modeling. This example describes the scenario when the initial basic reproduction number R 0 = 3.5 was reduced to R 0 = 1.1 by the quarantine measures at t 1 = 20 days. The number density of quarantine individuals continued to grow but at much lesser rate and acquired a local peak at t ≈ 48 days. Then the quarantine was terminated at t 2 = 60 days. The number densities of the infected and quarantine individuals immediately started to grow again and reached the new peaks that were much larger than those during the quarantine (row 5 in Table I ). As compared to the "no quarantine" scenario (row 3 in Table I ), the absolute heights of the peaks decreased, but the total number densities of quarantine individuals and fatal cases remained nearly the same. Thus, the quarantine scenario allows one to win time but does not seriously affect the total number of ill and deceased people, in the case when the fatality rate remains to be constant.

Our results are in line with the results of modeling presented in Ref. [22] . The greater the reduction in transmission, the longer and flatter is the epidemic curve, with the risk of resurgence when interventions are lifted to mitigate economic impact. The similar results were obtained when modeling the COVID-19 quarantine scenario for the Wuhan city, with a stochastic SEIR model fitted to the available statistical data [44] . The initial R 0 value equal to 2.35 (the median estimate) before the quarantine restrictions took place declined to R 0 ≈ 1.05 after the start of the quarantine.

Finally, we demonstrate the effect of nonzero resource activation parameter E. Figure 7 demonstrates the effect of the resource function for E = 0.1, with all other parameters listed in Table III . The nonzero resource activation parameter E has a profound effect on the number of fatal cases. This example illustrates the scenario when the economic subsystem has no resource to support the medical infrastructure, so that seriously ill (quarantine) individuals could not get the necessary medical services to overcome the infection. However, in contrast to the toy example considered in Section II, dynamic system (22) does not collapse, the resource function exhibiting a clearly identified minimum with a subsequent surge.

We proposed a novel approach for modeling the spread of epidemics. The spreading process within a selected social group is described by the dynamic equations which explicitly account for the dependence of characteristic reaction (infection) rates on the local population density in various social zones. The indirect channel of infection transmission via an intermediate environment or the so-called fomites was also taken into account. Moreover, the model accounts for the negative feedback in the "epidemic-economic resource" system. To this end, the equations describing the spreading process were supplemented with equations for the formal economic resource. The epidemic spread and the use of quarantine measures was demonstrated to be connected with economic losses, which in turn can lead to the collapse of the system, at least for a number of industries and/or social groups.

Appendix A: Model parameters (extended) Table II gives transmission frequencies ω j and rates β calculated by formulas (12) and (14) for various sets of social control parameters. Table III gives the full description of the model parameters and their estimates used in our computations.

In particular, the cloud-related parameters can be estimated as follows. The typical frequency Ω j of the pathogen transmission from the cloud j to the susceptible individual can roughly be estimated as

where Γ js is the average rate the susceptible individual contacts the cloud j, υ c is a typical volume of the pathogen transferred from the cloud to the susceptible individual per one contact, ∆ is some characteristic weight of one pathogen specimen that can be interpreted as the minimum portion ("quant") of the pathogen that can be transferred per one contact), and θ is the probability the transmission of this quant results in infection. The smaller the pathogen, the smaller is the quant and the more intensive is transmission (Ω j is higher). The parameters υ c , ∆, and θ are assumed to be independent of the particular cloud. The pathogen shedding rate σ j (infected-to-cloud) can roughly be estimated as

where Γ ij is the average rate the infected individual contacts the cloud j (number of coughs, sneezes, touches, etc. per unit time), n is the typical number of the pathogen quants ∆ transferred by the infected individual to the cloud per one contact, and V j is the total volume (capacity) of the cloud j. The parameter n is assumed to be independent of the particular cloud. Then the scaled indirect transmission frequencies ν j can be estimated as Table II βj Indirect transmission rate (infected-cloud-susceptible) in cloud j see Table IV βp Aggregate indirect transmission rate (infected-cloud-susceptible) 0. Tj Average time spent in location j see Table II ℓ0 Minimum possible distance between two individuals 1 m ℓc Correlation radius (the maximum transmission distance) 4 m ℓj Average distance between individuals in location j see Table II τis Characteristic time of becoming infected at close contact T /3 τp Average pathogen decay time outside the host 2 T 0. i0 Initial number density of infected individuals 10 −5 † for COVID-19 ‡ for cholera outbreak § for pandemic influenza * depends on medium, ambient temperature, and surface type * * 80% of COVID-19 cases are mild or asymptomatic

The dimensionless parameter

defines the transmission efficiency from the infected individual to the susceptible individual through the cloud j. Although the abstract parameter ∆ was eliminated by scaling (21) , the expression for χ j still contains the parameters that are hard to estimate from some physical principles. Therefore, this parameter can rather be estimated by fitting the model to some real statistical data.

In practice, it is selected by assuming that the indirect and direct routes transmissions have approximately the same likelihood, i.e. ν j ≈ ω j [32] . Table IV gives indirect transmission frequencies ν j and rates β j calculated by formulas (A3) and (23) for a par- Li et al. [49] Wuhan 2.5-2. 9 Wu et al. [56] Wuhan ticular set of social control parameters. The contact infected-to-cloud rates Γ ij , pathogen decay rates γ j , and transmission efficiencies χ j are assumed to be the same for each cloud and listed in Table III . The contact cloudto-susceptible rate Γ js is assumed to be inversely proportional to the squared social distance ℓ j between the individuals [as in Eq. (12) ]. Table V lists some estimates of the basic reproduction number R 0 derived from the statistical data on COVID-19 (literature data).

Epidemic processes in complex networks

Predicting the speed of epidemics spreading in networks

Small inter-event times govern epidemic spreading on networks

Rare and extreme events: the case of COVID-19 pandemic

Tail risk of contagious diseases

Macro-modelling and prediction of epidemic spread at community level

The mathematics of infectious diseases

SIS and SIR epidemic models under virtual dispersal

A Short History of Mathematical Population Dynamics

Mathematical Models in Epidemeology

Outbreaks in susceptibleinfected-removed epidemics with multiple seeds

Severe population collapses and species extinctions in multihost epidemic dynamics

English translation entitled "An attempt at a new analysis of the mortality caused by smallpox and of the advantages of inoculation to prevent it

Sur l'application du calcul des probabilitésà l'inoculation de la petite vérole

Daniel Bernoulli's epidemiological model revisited

The Prevention of Malaria

A contribution to the mathematical theory of epidemics

Contributions to the mathematical theory of epidemics. II. The problem of endemicity

Contributions to the mathematical theory of epidemics. III. Further studies of the problem of endemicity

Impact of the infectious period on epidemics

Reaction-diffusion processes and metapopulation models in heterogeneous networks

The spread and control of rumors in a multilingual environment

How will country-based mitigation measures influence the course of the COVID-19 epidemic?

Vaccination and the theory of games

The scaling of contact rates with population density for the infectious disease models

Merging economics and epidemiology to improve the prediction and management of infectious disease

Chemical Kinetics

Arrhenius Equation and Non-Equilibrium Kinetics: 100 Years Arrhenius Equation

The Theory of Rate Processes: The Kinetics of Chemical Reactions, Viscosity, Diffusion and Electrochemical Phenomena

Energetic processes in macroscopic fractal structures

Diffusion and Heat Transfer in Chemical Kinetics

On the theory of uniform flame propagation

Multiple transmission pathways and disease dynamics in a waterborne pathogen model

Global dynamics of cholera models with differential infectivity

Koopman, Dynamics and control of infections transmitted from person to person through the environment

Fomite-mediated transmission as a sufficient pathway: a comparative analysis across three viral pathogens

The Devil's staircase

Why is it difficult to accurately predict the COVID-19 epidemic?

Phase-adjusted estimation of the number of coronavirus disease 2019 cases in Wuhan

Mathematical modelling of the transmission dynamics of Ebola virus

Recurrent outbreaks of childhood diseases revisited: the impact of isolation

Identifiability and estimation of multiple transmission pathways in cholera and waterborne disease

A new epidemic model with indirect transmission

Epidemic models with heterogeneous mixing and indirect transmission

Early dynamics of transmission and control of COVID-19: a mathematical modelling study

Quantifying the routes of transmission for pandemic influenzas

Pathogenicity and transmissibility of 2019-nCoV A quick overview and comparison with other emerging viruses

Stability of SARS-CoV-2 in different environmental conditions

The incubation period of coronavirus disease 2019 (COVID-19) from publicly reported confirmed cases: estimation and application

Early transmission dynamics in Wuhan, China, of novel coronavirusвҐҮinfected pneumonia

Symptom progression of COVID-19

Response Team

Estimating clinical severity of COVID-19 from the transmission dynamics in Wuhan, China

Estimates of the severity of coronavirus disease 2019: a model-based analysis

Transmissibility of 2019-nCoV

Novel coronavirus 2019-nCoV: early estimation of epidemiological parameters and epidemic predictions

Preliminary estimation of the basic reproduction number of novel coronavirus (2019-nCoV) in China, from 2019 to 2020: A data-driven analysis in the early phase of the outbreak

Nowcasting and forecasting the potential domestic and international spread of the 2019-nCoV outbreak originating in Wuhan, China: a modelling study