key: cord-0904481-rll24klp authors: Gandzha, I.S.; Kliushnichenko, O.V.; Lukyanets, S.P. title: Modeling and controlling the spread of epidemic with various social and economic scenarios date: 2021-06-03 journal: Chaos Solitons Fractals DOI: 10.1016/j.chaos.2021.111046 sha: 6ddb2442402eea1be4010f5d8c5ffa3b2104d593 doc_id: 904481 cord_uid: rll24klp We propose a dynamical model for describing the spread of epidemics. This model is an extension of the SIQR (susceptible-infected-quarantined-recovered) and SIRP (susceptible-infected-recovered-pathogen) models used earlier to describe various scenarios of epidemic spreading. As compared to the basic SIR model, our model takes into account two possible routes of contagion transmission: direct from the infected compartment to the susceptible compartment and indirect via some intermediate medium or fomites. Transmission rates are estimated in terms of average distances between the individuals in selected social environments and characteristic time spans for which the individuals stay in each of these environments. We also introduce a collective economic resource associated with the average amount of money or income per individual to describe the socioeconomic interplay between the spreading process and the resource available to infected individuals. The epidemic-resource coupling is supposed to be of activation type, with the recovery rate governed by the Arrhenius-like law. Our model brings an advantage of building various control strategies to mitigate the effect of epidemic and can be applied, in particular, to modeling the spread of COVID-19. The spread of contagions, diseases, information, rumors, ideas, or concepts share many similarities. In most cases, such spreading processes are governed by similar models [1] [2] [3] [4] [5] [6] , which can vary in terms of complexity or approximations. The simplest dynamical models (see, e.g., Refs. [7] [8] [9] [10] [11] [12] [13] [14] ) have their origin in classical works [15] [16] [17] [18] [19] . More complex models include stochastic effects [20] [21] [22] , spatial flows (e.g., diffusion [23, 24] ), or allow for nontrivial spatial structure or topology [25] [26] [27] [28] [29] [30] . One common feature in the majority of these models is the presence of kinetic coefficients or parameters that characterize the probability of elementary processes (reactions) per unit time. For instance, contagion, infection, or disease spreading is associated with scattering of infected individuals on the spatiotemporal fluctuations of population density. The probability of being infected is not a monotonous function of population density [31, 32] . The probability of scattering of infected individuals on susceptible individuals is lower at low densities, which results in a smaller total scattering cross-section. At higher densities, the mobility of spreaders decreases suppressing 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., an epidemic. The kinetic description of spreading is generally quite a complicated problem that should take account of the internal state of spreaders, which changes in the course of collisions, the presence of spatial inhomogeneities, and non-equilibrium properties of the system. The corresponding kinetic coefficients essentially depend on the local population density and characteristic time scales of collision processes. The dynamics of any spreading process such as an epidemic is determined by individual peculiarities of people in a selected social group, by the type and mechanism of infection, etc. It can also be affected by the ways the process is controlled and influenced. A dramatic example is the spread of COVID-19, when different countries and governments resorted to different control strategies and quarantine measures [33] [34] [35] . Such measures can be controlled by the choice of certain parameters accessible to society (e.g., working day duration, average density of people in public places, frequency of disinfection, intensity of transport communications, etc.). The problem of strategy selection reduces to problems of optimal control theory for feedback systems [36] [37] [38] or to the theory of games in a more general case [39] . Spreading processes are subject to instabilities. On the one hand, the instability with respect to initial conditions and local fluctuations can produce a critical number of secondary cases and trigger a collective chain reaction like an epidemic. On the other hand, the parametric instability of spreading processes is often attributed to certain values of kinetic coefficients. However, 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 like population density, characteristic time scales, etc. For example, the recovery rate is generally determined by the quality of provision with medical services and food, apart from the individual peculiarities of the given member of population (see, e.g., Refs. [40] [41] [42] [43] [44] [45] ). The quickest recovery depends on the cost of medical services and the bare subsistence level of consumption E, as well as on the availability of collective economic resource ρ. Since the cost of services is fixed, the service is terminated if there is no sufficient resource ( ρ E). In other words, the parameter E serves as the height of some energy barrier (activation energy) peculiar to the given system. As was shown in Ref. [46] , the recovery rate can have an activation-type (Arrhenius-like) dependence, ∝ exp (−E/ρ) , similar to the temperature dependence of common activation processes with activation energy E [47, 48] . In this case, the economic resource ρ formally plays the role of effective market temperature and the minimum level of resource consumption is associated with activation energy E. This observation follows from the fact the effective market temperature can be associated with the average amount of money or average income per individual (or economic agent) and the equilibrium distribution of income or money for single economic agents is governed by the exponential Boltzmann law (at least for low and middle income classes) [49] [50] [51] [52] . The activation process implies that the system can exhibit the so-called explosive (or catastrophic) instability [53] . 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 ultimately to a thermal explosion, which is described in the framework of the Zel'dovich-Frank-Kamenetskii theory [54] [55] [56] . The fight against the epidemic involves similar catastrophic processes. The spread of epidemic and the associated quarantine measures result in the reduction of the collective resource ρ. When the resource is depleted, the quality of medical services drops and the recovery rate goes down. As a result, the number of active members in population decreases. This, in turn, leads to a further reduction of the collective resource, with the level of income needed for the basic survival being lower and lower. Such a scenario can finally result in the ultimate collapse of the system-the effect opposite to thermal explosion [46] . In this paper, we address some of the points mentioned above. In other words, we are interested not only in modeling the spread of epidemic as an example of a 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 a quite simple dynamical feedback model. This model has the following two important features: (I) To take into account the fact that an individual can stay in different environments or public places (e.g., transport, shop, work) characterized by different densities of surrounding individuals, it is convenient to consider some social group and its average daily cycle ( T ). For instance, we can introduce the average daily time spans for which the individual stays 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 in each of the locations j. We also make a qualitative assessment of the dependence of the transmission rate and other kinetic coefficients on the density of individuals or average distance between them in order to clearly identify the set of control parameters { T j , j } and their effect on the spreading process. Finally, we take into consideration a possibility of indirect transmission through the medium (e.g., contaminated water [57, 58] ) or via the so-called fomites [59, 60] . The indirect transmission rate depends on the same set of control (II) The interplay between the epidemic and the economic resource available for the selected social group is taken into account by means of the activation-type dependence of characteristic recovery rates. To this end, we use the simplest form of the equation for the dynamics of the collective resource ρ (such as in Ref. [46] ). Note that for the sake of simplicity we do not take into account the contribution of local density fluctuations but consider the averaged rate constants assuming that the flows of individuals between different spatial locations are counterbalanced during the selected daily cycle. In other words, we do not consider the phase transitions induced by noise or fluctuations (see, e.g., Ref. [61] ). The paper is organized as follows. In Section 2 , we estimate the direct and indirect transmission rates of contagion in terms of basic social control parameters such as the average distance between the individuals and the average time spent in a certain location. In Section 3 , we consider an extended suspected-infectedquarantine-recovered (SIR/SIQR) model that takes into consideration the above-discussed points (I) and (II). In Section 4 , we provide several examples of numerical modeling. In particular, we demonstrate the effect of indirect transmission, the effect of quarantine measures, and the effect of limited resource. Conclusions are drawn in Section 5 . As mentioned above, the spread of contagion such as infection can be controlled only by the parameters accessible to society. For example, the average distance j between the individuals and the average time T j spent in a certain location j can be such control parameters. Here we estimate the infection transmission rates in terms of these parameters in the case of direct transmission (infected-to-susceptible) and indirect transmission (through the environment). It can naturally be assumed that the probability of being infected, which generally depends on the distance between the infected and susceptible individuals, is described by a damped function with a characteristic correlation radius c (decay length). For simplicity, we suppose that the susceptible individual and the infected individual can interact only when they fall into the vicinity of radius c with each other and that the probability of elementary transmission per unit time in this case does not depend on distance and equals ω. We also assume that for any susceptible individual the transmissions from two different inf ected individuals are independent events. Note that in our simplified consideration the correlation radius c can be associated with a socially safe distance. Next we suppose that for each social location j the average population density c j is not affected by any social conditions and is solely determined by the average distance between individuals, The average number of individuals falling into the Table 1 Social control parameters and transmission rates for 0 = 1 m, c = 4 m, and τ is = T / 3 ( T = 1 day). circle of radius c around the susceptible individual in each of the locations j is N j ≈ c j π 2 c . Then the probability for the susceptible individual to be infected over time interval t can be estimated as Here κ ≈ i = I/N is the probability that the given individual in the vicinity of the susceptible individual is infected, N being the total number of individuals, I being the total number of infected individuals, and i being the number density of infected individuals. Summing over all the susceptible individuals S, we obtain an increment of the number of infected individuals, I = S P j . Dividing by the total population N and the time interval t, we get a characteristic growth rate of the number density of infected individuals: where s = S/N is the number density of susceptible individuals. Thus, the direct (infected-to-susceptible) transmission rate of infection in the selected area with population density c j is estimated as The probability ω of elementary transmission per unit time can be estimated in terms of the minimum possible distance 0 between two individuals (the so-called close contact distance). Assume that the susceptible individual is infected with probability equal to unity if he stays at the distance 0 from the infected individual for some characteristic time τ is . Then the probability of becoming infected at the minimum distance per unit time is propor- is . The parameter a can be estimated as the ratio of the contact area between the susceptible individual and the infected individual to the total area determined by the correlation radius c : a ≈ ( 0 / c ) 2 . Thus, we have Let us assume that each individual can stay in N ω different locations during some typical period 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 the following formula: where T j are the characteristic daily time spans spent in each of the locations j, with As an example, consider four different social locations ( N ω = 4 ). Location 1 refers to staying at home (limited social contacts), 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. Table 1 gives the transmission rates ω j and β calculated for various sets of social control parameters T j and j . 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. It can be seen that our rough estimates given by formulas (2) and (3) 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 ( ≈ 0 . 6 T −1 ) derived from the statistical data for the COVID-19 epidemic in the Wuhan city [62] . Note that the ω j estimate given by formula (2) is very rough. In fact, the transmission rate depends on fluctuations in the number of individuals falling into the circle of radius c . The maximum number N m of individuals that can fall into the circle of radius c is determined by the closest packing of the region r < c by hard spheres with radius equal to the the minimum possible distance 0 between two individuals. For example, considering a triangular lattice with lattice constant 0 , we have Let ϑ be the probability that the lattice node is occupied. Then the average number of individuals falling into the region r < c is N = ϑN m , and its mean-square fluctuation is estimated as Finally, the relative error introduced by such fluctuations in the ω j estimate can therefore be estimated as It should also be noted that the characteristic scattering crosssection given by relation (1) can have a more complex dependence on density, ∼ c α ( 1 ≤ α ≤ 3 ), and be governed by a nonlinear incidence function other than is [31, 32, [63] [64] [65] [66] . In addition to the direct transmission of infection by direct contact of susceptible individuals with infected individuals, the infection caused by a pathogen (e.g., virus) can be transmitted indirectly either through the medium (e.g., contaminated water [57, 58] ) or via intermediate objects (like hands, counters, doorknobs, etc.) often called fomites [59, 60] . Such an intermediate medium or fomites will further be referred to as "cloud". We associate a separate cloud with each social location j. Let d j be the pathogen density per individual in the cloud j. Then a characteristic growth rate of the number density of infected individuals attributed to the indirect transmission of the infection from the cloud j to susceptible individuals can be calculated in the same way as in the case of direct transmission, namely Here the parameter j describes the typical (indirect) rate of pathogen transmission from the cloud j to a susceptible individual. The estimates of indirect transmission rates for various sets of control parameters T j and j are provided in Appendix A . We subdivide the selected social group into five compartments: susceptible (S), infected (I), quarantined (Q), recovered (R), and deceased (F). Susceptible individuals (S) are those who are at risk to be infected. 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 disease (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 quarantined individuals (Q). The quarantined 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 operator ∂ t stands for the derivative with respect to time t. The functions s (t) , i (t) , q (t) , and r(t ) describe the number densities of the susceptible, infected, quarantined, 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 2 ). The direct transmission rates ω j = ω j ( j ) in each of the locations j are defined by formula (2) , and the indirect transmission rates The characteristic daily time spans T j in each of the locations j are given in Table 1 , the total time spent in all the locations being constant each day and equal to the day duration T [see formula (4) ]. Note that in general each location j should also be characterized by a certain local temperature that determines the mobility of spreaders (see, e.g., Ref. [31] ). However, we do not consider such an additional variable in this work for the sake of simplicity. Each of the functions d j (t) describes the pathogen density 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: (6e) The first term, σ j i , describes the pathogen shedding by the infected individuals into the cloud. The coefficient σ j = σ j ( j ) defines the average pathogen shedding rate in the cloud j. The second term, γ j d j , describes the pathogen decay in the cloud due to natural inactivation, decontamination, or other routes. As a matter of fact, Eq. (6e) represents an equation for contagion in the location j. The function ρ(t) represents the average resource associated with the average amount of money or income per individual in the selected social group. The resource balance equation is written as follows [46] The acquisition of this resource per unit time is proportional to the number density of working (active) individuals (the quarantined individuals are assumed to be not working). The parameter 4 ) formalizes the resource amount acquired by them per unit time. It is proportional to the average working time T 4 for the given social group and the corresponding average population density c 4 ∝ ( 4 ) −2 . In line with Ref. [67] , we refer to this parameter as acquisition rate. 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 parameter represents a resource source (constant resource inflow into the system from some external reservoir) or a resource sink (constant resource outflow from the system). When > 0 , resource is fed into the system (e.g., in the form of subsidies) from some external source, e.g., a central bank or central government. When < 0 , resource flows out from the system, e.g., in the form of infrastructure expenses, depreciation, rent, interest payments, or other fixed expenses. The recovery process is governed by the general economic situation characterized by a certain minimum level of resource consumption E, which reflects the cost of medical and other essential life services. It defines the minimum amount of the consumed resource needed for the recovery of individuals in the selected social group, therefore implying the existence of some "energy" barrier for their recovery. The recovery rates, i.e., the coefficients at the function q (t) in Eqs. (6a) and (6d) , are supposed to have an activation-type (Arrhenius-like) dependence, ∝ exp (−E/ρ) , similar to the temperature dependence of common activation processes with activation energy E. The number density of fatal cases is described by the function f (t) that is determined by the equation The above equations are supplemented with the following initial conditions i 0 being the initial number density of infected individuals and ρ 0 being the initial resource value. The total number density of individuals is assumed to be constant: Effectively, our model is a combination or extension of SIR-like (susceptible-infected-recovered) epidemic models such as SIQR (susceptible-infected-quarantined-recovered) [68] , SIWR/SIVR/SIRP (susceptible-infected-recovered-pathogen) [57, [69] [70] [71] , SIRD (susceptible-infected-recovered-deceased) [72, 73] , and EITS (environmental infection transmission system) [59] . The extension to multiple groups is given in Appendix B . 1. The population subsystem is closed (no migration outside the selected population group and no one is added to the population). 2. Natural demography is ignored. 3. A uniform spatial distribution of people is assumed. The pathogen distribution in each cloud is assumed to be uniform. ω j , j , is , ir , iq , qs , qr , and q f between the compartment levels are given in Table A.1 . 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-exposedinfected-recovered) model was demonstrated to have no practical advantage as compared to the SIR model [62] . 5. Quarantined individuals do not infect others. 6. There is no pre-existing immunity in susceptible individuals. 7. There is no loss of immunity by recovered individuals. The full description and indicative values of the parameters of Eqs. (6) are given in Table A .1 (see Appendix A ). Figure 1 shows a schematic diagram that depicts transitions between different compartments and identifies the corresponding transition rates. In particular, the parameter is describes a rate 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 a 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 A .1 ). The rate constant i is defined as follows: where τ i is the characteristic pathogen incubation period. The parameters qr and qs describe the rates at which the quarantined individuals recover with or without the acquired immunity. The parameter q f describes the fatality rate in the case of unlimited resource ( E ρ). In this case, the quarantined 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 ρ (no resource to fight against the epidemic), the fatality rate is at its maximum and is equal to the rate constant q = qr + qs + q f = τ −1 q for the quarantined individuals. Here τ q is the mean quarantine/hospitalization period. This case refers to the full collapse of the medical system when the quarantined individuals get no medical help or treatment. Each of the parameters j describes the typical rate of pathogen transmission from the cloud j to susceptible individuals. 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 . Note that the instantaneous number density of quarantined (ill) individuals q (t) is not often a convenient indicator for practical applications. The integral number density of quarantined individuals for a certain period of time can be used instead. It is calculated as In the case of unlimited resource ( E = 0 ), the integral number density of quarantine individuals in the end of epidemic ( t → ∞ ) is proportional to the total number density of fatal cases, where η is the probability of the fatal scenario for a quarantined individual. For our further analysis, we first find the stationary solution to Eq. (6e) 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. (6) can be rewritten as Here β is the direct transmission rate (infected-to-susceptible) given by formula (3) . The parameters are defined in terms of the scaled indirect transmission rates (infected-cloud-susceptible) via each of the clouds j. Eqs. (12) supplemented with Eqs. (6f) and (6g) possess two equilibrium points. The first one is the disease-free equilibrium The second one is the endemic equilibrium where the parameter is the basic reproduction number and is the integral indirect pathogen transmission rate via all the clouds. In general, the basic reproduction number R 0 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 state of endemic equilibrium. The first two equations of system (12) 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 → ∞ : ln Eq. (19) can be used to control the accuracy of the numerical integration of system (12) . In the examples considered in the next Section, we used the fourth-order Runge-Kutta method to integrate Eqs. (12) numerically with step t = T / 10 sufficient to achieve the reasonable accuracy. In the case of μ = 1 (no loss of immunity), our numerical estimate of s ∞ coincided with the value given by Eq. (19) to an accuracy of 10 −10 . Now we consider particular examples to demonstrate various effects described by our model. All the control parameters { T j , j } are selected according to Table 1 for the selected social scenario (casual, soft quarantine, strict quarantine). All other parameters are either estimated in terms of control parameters or listed in Table A. 1 . In particular, the integral direct and indirect transmission rates β and β p are calculated in terms of parameters { T j , j } by formulas (3) and (18) , respectively. First we focus on the case when there is no resource depletion ( ρ E), so that the resource activation barrier could be ignored ( E = 0 ). We illustrate the "patient zero" phenomenon, demonstrate the effect of indirect transmission, and model a quarantine scenario. Next we consider an example of a social group with limited resource ( E = 0 ). Figure 2 a shows the number densities of susceptible, infected, quarantined, 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 ). The case i 0 = 10 −7 corresponds to an initial density of one per 10 million, and the case i 0 = 10 −5 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 2 ). 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. 3 ) . Figure 2 shows the number densities of susceptible, infected, quarantined, and recovered individuals in the cases when there is no indirect transmission [panel (a)] and when there is such a transmission [panels (b,c) ]. The model parameters were selected according to Table A .1 (see Appendix A ). 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 might become larger and shift to shorter times ( Table 2 ). 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 larger are i max , q max , and (q ) ∞ . Thus, to reduce the epidemic peak and to slow the epidemic down, one should reduce R 0 . According to formula (17) , this can be achieved by reducing the transmission rates ω j and, therefore, the integral transmission rate β. Such measures are usually referred to as quarantine. To model the quarantine scenario, we assumed the following form of the transmission rates: where t 1 is the moment when the quarantine starts, t 2 is the moment when the quarantine ends, β(0) and β p (0) are the direct and indirect transmission rates during the period when there is no quarantine, and β and β p are the transmission rates during the quarantine period. Figure 4 demonstrates an example of the quarantine scenario when the initial basic reproduction number R 0 = 3 . 5 was reduced to R 0 = 1 . 1 by quarantine measures at the moment t 1 = 20 days. In particular, this can be achieved by increasing the average distance j between the individuals in transport and other social locations and by reducing the average time T j spent in these locations (see Table 1 in Sect. 2). The number density of quarantine individuals continued to grow during the quarantine period but at much lesser rate and acquired a local peak at t ≈ 48 days. Then the quarantine was terminated at the moment t 2 = 60 days. The number densities of the infected and quarantined individuals immediately started to grow again and reached the new peaks that were much larger than those during the quarantine (row 5 in Table 2 ). As compared to the "no quarantine" scenario ( Fig. 2 b, row 3 in Table 2 ), the absolute heights of the peaks decreased, but the integral number densities of quarantine individuals q and, therefore, 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 by the end of epidemic, in the case when the mortality rate remains to be constant. The second (post-quarantine) peak in the number densities of infected and quarantined individuals clearly illustrates the effect known as the second wave of the epidemic, which was in particular observed in the case of COVID-19 epidemic in many countries [74, 75] . More complex dynamics and intervention measures can result in successive epidemic waves [76, 77] . Our results are in line with the results of modeling presented in Ref. [33] . 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 [78] . The pre-quarantine R 0 value equal to 2.35 (the median estimate) dropped to R 0 ≈ 1 . 05 after the start of the quarantine. Here we demonstrate the effect of nonzero resource activation parameter E. Let us rewrite resource Eq. (6f) in terms of dimen- where ρ (0) = (G (0) + ) / ρ is the equilibrium resource value when there is no epidemic (so that (0) ≡ 1 ), s = /G (0) is the number density of active individuals that would have to be working to acquire the resource amount equal to | | (per unit time), G (0) is the resource acquisition rate when there is no epidemic and the system is in equilibrium, and parameter is the time-dependent coefficient that allows for variations in the resource acquisition rate during the spread of the epidemic. For the sake of simplicity, we will limit our consideration to the case s = 0 . Resource Eq. (21) is coupled to Eq. (12) , where the exponential factor exp (−E/ρ) needs to be identically rewritten as exp (−E/ ) . The parameter E ≡ E/ρ (0) is the resource activation level E normalized to the equilibrium resource value. Nonzero E reduces the recovery rates of quarantined individuals, qr exp (−E/ ) and qs exp (−E/ ) , so that there would be more fatal cases as compared to the case of unlimited resource ( E = 0 ). Figure 5 demonstrates the effect of limited resource in the case of E = 0 . 1 , with all other parameters selected such as in the example shown in Fig. 2 b . Nonzero E has a profound effect on the number of fatal cases, which is nearly 6 times larger than in the case of unlimited resource ( Fig. 5 a) . Such a manyfold increase in the number of fatal cases is caused only by a 5% drop of resource ( Fig. 5 b) . Resource began to decline as the number of quarantined individuals grew up ( Fig. 5 a) , since the quarantined individuals are supposed to be passive and not acquiring the resource [see Eq. (21) ]. After the number of quarantined individuals passed through its peak, resource passed through its minimum and started to increase towards its initial value. This example illustrates a scenario when the economic subsystem has a limited capacity to support the medical infrastructure, so that seriously ill (quarantined) individuals could not get the necessary medical help to overcome the infection. Table 1 ) with t 1 = 20 T and t 2 = 60 T are β(0 The resource acquisition rate G was assumed to be constant ( k ≡ 1 ) in the above example. The parameter G is determined by the average number of working hours per working individual and by the working efficiency. The average number of working hours is proportional to the control parameter T 4 , which can be reduced in the case of quarantine, as discussed in Section 2 (see Table 1 ). Figure 6 demonstrates the effect of two different quarantine scenarios on the number density of fatal cases. When there is no quarantine, the number density of fatal cases is the same as in Fig. 5 a for the case of E = 0 . 1 . In the case of soft quarantine, the transmission rate β goes down ( Table 1 ) , the basic reproduction number R 0 becoming smaller and the epidemic spreading being less intensive. The soft quarantine scenario does not affect the average number of working hours per working individual, and the resource acquisition rate G remains the same as in the case of no quarantine ( k = 1 , here we ignore the minor effect of control parameter 4 ). The number of fatal cases grows much slower during the quarantine, but it gradually goes back towards its value in the case of no quarantine after the quarantine is terminated. In the case of strict quarantine, the transmission rate β is yet smaller ( Table 1 ) , mainly because this quarantine scenario also affects the total number of working hours, which is twice as small as compared to the case of soft quarantine. As a result, the resource acquisition rate goes down ( k ≈ 0 . 43 , here we also take into account the effect of control parameter 4 , see Eq. (22) ) during the quarantine. In the short run, the number of fatal cases in the case of such strict quarantine measures substantially declines as compared to the cases of soft quarantine and no quarantine ( Fig. 6 ) . However, it starts to rapidly increase after the quarantine is terminated and eventually becomes larger than it was in the case of no quarantine. This example demonstrates that strict quarantine measures that affect the general economic situation may have serious negative social outcomes in the case of limited economic resource. We considered a dynamical model for describing the spread of epidemics. The spreading process within a selected social group is governed by the equations that explicitly take into account the dependence of characteristic transmission rates on the local population density in various social zones. The indirect channel of transmission via an intermediate environment or the so-called fomites was also taken into consideration. A negative feedback between the infected population size and a collective economic resource associated with the average amount of money or income per individual was introduced to describe the socioeconomic interplay. The epidemic spread and the use of quarantine measures was demonstrated to be connected with economic losses, which in turn could aggravate the negative outcomes of the epidemic. Note that in this paper we did not take into consideration the local spatial dispersion inside the selected spatial zones and any transitions from one zone to another (network models). As a matter of fact, averaging over the temporal period (daily cycle) of a quasi-periodic system allowed us to roughly eliminate fast variables and to avoid the need of accounting for the probabilities of transitions between spatial zones. The model presented here can be used to model the COVID-19 epidemic for particular social groups and regions. However, matching to real statistical data was not the direct purpose of this work and can be performed as the further development of this study. The model can also be applied to describing other spreading processes, such as the spread of information, rumors, ideas, or concepts. The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. O.K. was partially supported by a grant for research groups of young scientists from the National Academy of Science of Ukraine (Project No. 0120U100155 ). We thank Prof. B.I. Lev for fruitful discussions. Table 1 in Section 2 gives the direct transmission rates ω j and β calculated by formulas (2) and (3) for various sets of social control parameters. Table A.1 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 (indirect) rate j of pathogen transmission from the cloud j to a 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 the transmission ( j is higher). The parameters υ c , , and θ are assumed to be independent of the particular cloud. The rate constant js is assumed to be proportional to the average population density c j ∝ ( j ) −2 in the corresponding cloud. The pathogen shedding rate σ j (infected-to-cloud) can roughly be estimated as where i j 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. The rate constant i j is assumed to be proportional to the average population density c j in the corresponding cloud. Then the scaled indirect transmission rates ν j can be estimated The dimensionless parameter defines the transmission efficiency from infected individuals to susceptible individuals through the cloud j. Although the parameter was eliminated by scaling (11) , 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 of transmissions have approximately the same likelihood, i.e. ν j ≈ ω j [57] . Resource inflow or outflow per unit time 0 i 0 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 Scaled indirect (infected-cloud-susceptible) transmission rates ν j and β j in the case of the casual epidemic scenario (see Table 1 ). The aggregate indirect transmission rate βp is given by formula (18) . T j , h j , m The pathogen decay rates γ j and transmission efficiencies χ j are assumed to be the same for each cloud and listed in Similarly to the basic SIR model [7, 9] , our model can easily be extended to the multigroup formulation, e.g., with subdivision by age. The dynamics of each group n is governed by the following equations: The pathogen dynamics in the cloud j is given by an equation Analysis and control of epidemics. A survey of spreading processes on complex networks Critical behaviors in contagion dynamics Efficient spread-size approximation of opinion spreading in general social networks Predicting the speed of epidemics spreading in networks The spread and control of rumors in a multilingual environment Macro-modelling and prediction of epidemic spread at community level The mathematics of infectious diseases SIS and SIR epidemic models under virtual dispersal Mathematical modeling of COVID-19 transmission dynamics with a case study of Wuhan The first 100 days: modeling the evolution of the COVID-19 pandemic A short history of mathematical population dynamics An introduction to mathematical epidemeology Mathematical models in epidemeology An attempt at a new analysis of the mortality caused by smallpox and of the advantages of inoculation to prevent it" in: Bradley L, Smallpox Inoculation: An Eighteenth Century Mathematical Controversy (Adult Education Department Heesterbeek JAP 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 SBDiEM: a new mathematical model of infectious disease dynamics The dynamics of entropy in the COVID-19 outbreaks Reaction-diffusion processes and metapopulation models in heterogeneous networks Localized outbreaks in an S-I-R model with diffusion Epidemic processes in complex networks Outbreaks in susceptible-infected-removed epidemics with multiple seeds Small inter-event times govern epidemic spreading on networks Spreading of infections on random graphs: a percolation-type model for COVID-19 Modeling the spatiotemporal epidemic spreading of COVID-19 and the impact of mobility and social distancing interventions Epidemic growth and Griffiths effects on an emergent network of excited atoms The scaling of contact rates with population density for the infectious disease models Law of mass action and saturation in SIR model with application to coronavirus modelling How will country-based mitigation measures influence the course of the COVID-19 epidemic? Creating and applying SIR modified compartmental model for calculation of COVID-19 lockdown efficiency Predicting the dynamical behavior of COVID-19 epidemic and the effect of control strategies Merging economics and epidemiology to improve the prediction and management of infectious disease Optimal allocation of resources for suppressing epidemic spreading on networks Controlling percolation with limited resources Vaccination and the theory of games Dynamics of an SIR epidemic model with limited medical resources revisited Disease-induced resource constraints can trigger explosive epidemics Nontrivial resource requirement in the early stage for containment of epidemics Unravelling the myths of R 0 in controlling the dynamics of COVID-19 outbreak: a modelling perspective Dynamic analysis of a mathematical model with health care capacity for COVID-19 pandemic Epidemic spreading under infection-reduced-recovery Epidemic-driven collapse in a system with limited economic resource A toy model for the epidemicdriven collapse in a system with limited economic resource The theory of rate processes: the kinetics of chemical reactions, viscosity, diffusion and electrochemical phenomena Arrhenius equation and non-equilibrium kinetics: 100 years Arrhenius equation. Leipzig: BG Teubner Evidence for the exponential distribution of income in the USA Colloquium: statistical mechanics of money, wealth, and income Universal patterns of inequality Exponential structure of income inequality: evidence from 67 countries Explosive and nonexplosive onsets of instability The temperature distribution in a reaction vessel and the stationary theory of thermal explosions Energetic processes in macroscopic fractal structures Kinetic effects in thermal explosion with oscillating ambient conditions Multiple transmission pathways and disease dynamics in a waterborne pathogen model Global dynamics of cholera models with differential infectivity 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 Noise-induced transitions: theory and applications in physics, chemistry, and biology Why is it difficult to accurately predict the COVID-19 epidemic? Chemical reaction models for non-equilibrium phase transitions Influence of nonlinear incidence rates upon the behavior of SIRS epidemiological models Dynamical behavior of epidemiological models with nonlinear incidence rates A disease transmission model in a nonconstant population Interplay of sources of stochastic noise in a resource-based model 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 Analysis and forecast of COVID-19 spreading in China, Italy and France Characterization of the COVID-19 pandemic and the impact of uncertainties, mitigation strategies, and underreporting of cases in South Korea, Italy, and Brazil Second wave COVID-19 pandemics in Europe: a temporal playbook COVID-19: data-driven dynamics, statistical and distributed delay models, and observations Seoane-Sepúlveda JB . A SIR-type model describing the successive waves of COVID-19 Comparison of the second and third waves of the COVID-19 pandemic in south korea: Importance of early public health intervention Early dynamics of transmission and control of COVID-19: a mathematical modelling study Phase-adjusted estimation of the number of coronavirus disease 2019 cases in Wuhan Quantifying the routes of transmission for pandemic influenzas Stability of SARS-CoV-2 in different environmental conditions Pathogenicity and transmissibility of 2019-nCoV-a quick overview and comparison with other emerging viruses 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 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 Report 3 of the Imperial College London COVID-19 Response Team 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