On March 22, 2020 according to the official data presented at the Worldometer [17] for COVID-19, the number of infected individuals exceeded 10 thousand in seven countries. The respective growth rate of the number of infected is shown in Table 1 .
The first column represents the number of days needed to pass from 100 to 1000 infected, the second column from 1000 to 10000. If the number of days is the same in the two columns, then the growth is exponential, as it is predicted by conventional epidemiological models (SIR) for the beginning of epidemic when the number of uninfected (susceptible) individuals can be considered as approximately constant. In Italy, Iran, and France, the number of days in the second column is larger than in the first column. In Italy main restrictions were adopted beyond ten thousand initial infected. In Germany this difference is small at the moment, so it is early to make conclusions. In USA and Spain the number of days is the same and the growth rate is slightly larger in Spain. China  -7  Italy  6  10  USA  8  8  Spain  7  7  Germany  8  9  Iran  5 10 France 7 11 Table 1 . Tenfold time of the number of infected individuals in the country where their total number exceeds 10 thousand on March 22. The first column shows the number of days required to pass from 100 to 1000 infected, the second column from 1000 to 10 000 infected.
These data should be taken with caution because they can be strongly influenced by the number of effectuated tests. According to some media reports, there is a shortage of testing facility and/or shortage of testing kits in many countries, and the accurate number of infected individuals can be much larger [7] . Dynamics of the number of infected individuals can be influenced by various factors, in particular, by their spatial distribution and nature of social mixing. Appearance of new focuses of infection can essentially change the growth curves. For example, at the whole world level, development of the epidemic in Europe started when it was mostly at decaying mode in Asia. Instead of the standard growth-decay curves, this leads to growth-decaygrowth curves. Some elements of this more complex dynamics where the number of newly reported daily cases decreased during several days was also observed in some countries.
If we suppose that the distribution of infected individuals is uniform over space, then the SIR model predicts the exponential growth at the beginning of epidemic. Some time later, the growth rate decreases due to the decrease in the number of susceptible individuals, and the number of infected individuals decays at the end of epidemic. There are recent works with various modifications of the epidemiological models trying to take into account some specific features of the coronavirus pandemic [6, 9, 16] .
Analysis of the data on seasonal influenza shows that dynamics of epidemic spread can include the period of accelerated growth ( Figure 1 ). In the beginning of epidemic, the growth rate is exponential. At the next stage it continues to grow exponentially but with an increased growth rate. It can be considered as two different exponential functions where the second one has a larger exponent than the previous one. A possible explanation of this effect is related to the presence of sub-populations with different response on the infection. Simplifying the situation, we can consider two sub-populations, one of them with a weak immune response and another one with a strong immune response. The strength of immune response is determined by the rate of production of antigen specific immune cells. We will define it below. At the beginning of epidemic, infection propagates among the first subpopulation. When the epidemic reaches certain level, the second sub-population becomes also involved increasing the number of susceptible individuals.
Transition between the two stages of epidemic growth can be related to the number of infected individuals. If it is small enough, then a healthy individual rarely meets infected individuals, and the initial viral load remain small. People with a strong immune system can eliminate infection without developing disease symptoms, while people with a weak immune system can fall sick. At the later stages of epidemic, when the number of infected individuals is sufficiently large, an uninfected individual can cross several infected individuals during a short period of time. Therefore, the initial viral load becomes larger, and people with a strong immune response can also become ill.
In order to describe this two stage epidemic growth, we develop in this work an immunoepidemiological model combining the individual level and the population level. This model is generic but it can take into account specific features of a particular infection or the details of epidemiological situation (like quarantine, social distancing, self-isolation, etc.). In the next section, we will discuss some epidemiological data on influenza epidemic. Section 3 is devoted to a model of immune response where we show how people with weak and strong immune system can react on a given initial viral load. Immuno-epidemiological model is studied in Section 4.
There are large amount of data available for the seasonal influenza epidemic. A typical graph of the number of cases is shown in Figure 1 (left, blue curve) [5] . Growth rate of the number of positive cases between weeks 1 and 6 is exponential. The red curve shows the approximation of the blue curve with the function 20.5 exp(0.2x), where x is the number of weeks. The same curves are shown in Figure 1 (right) in the logarithmic scale. The red curve here is a straight line. The weeks 7-10 show different growth rate with some acceleration in growth. It also corresponds to an exponential growth but with a different exponent. Some other data are presented in Appendix 1. Figure 1 : Data on seasonal influenza in Europe in 2019-2020 [5] . The left graph shows the number of influenza-positive specimens from non-sentinel sources (blue line) and its approximation by an exponential (red line). The right graph shows the same curves in the logarithmic scale.
Conventional epidemiological models like SIR (Susceptible, Infected, Recovered, see, e.g., [11] ) and their numerous variants describe the density of infected individuals I by the typical equation
where S and I are the densities of susceptible and infected individuals, and λ is a positive constant. The first term in the right-hand side of this equation characterizes the appear-ance of new infected individuals due to their contact with susceptible. The second term corresponds to the decrease of I due to the recovery or death of infected individuals.
In the beginning of the infection spreading, the number of infected and recovered individuals are much less than the number of susceptible, so that we can approximate S by a constant S ≈ S 0 . Using this approximation, we obtain a linear differential equation with constant coefficients. Its solution is given by the function I(t) = I 0 exp(αt), where I 0 is the number of infected individuals at the initial moment of time, and α = λS 0 −ν. If α > 0, then I(t) exponentially growth. The same condition can be written as R 0 > 1, where R 0 = λS 0 /ν is the basic reproduction number.
Thus, according to equation (2.1), growth of the number of infected individuals is exponential in the beginning of epidemic, and it slows down later when S decreases. Exponential growth in the beginning corresponds to the data in Figure 1 . However, it does not describe accelerated growth during the second stage of epidemic. A possible explanation of this acceleration is that in the beginning the epidemic, disease spreads mainly among the people with weak immune response. At the second stage, people with strong immune response become also exposed. In order to investigate this hypothesis, we consider a model of immune response in the next section.
We consider a qualitative model of immune response suggested in [2] :
where v and c are the concentrations of viruses and immune cells within the organism. The first term in the right-hand side of equation ( The maximal level of viral concentration is determined by the initial viral load v 0 , the virulence of infection w, and the strength of immune response s. For a given virulence of infection w > 1, the maximal viral load increases with for larger initial viral loads and decreases for stronger immune response. If we set W = w ln w−w+1, then v * = v 0 exp(W/s).
Suppose that there are two levels of immune response s 1 and s 2 , with s 1 < s 2 . The first one corresponds to weak immune response and s 2 to strong immune response. For a given initial viral load v 0 , the corresponding maximal virus concentrations are v * 1 and v * 2 , where v * 1 > v * 2 . Disease development depends on virus concentration in a threshold way. If its level exceeds some critical value, which is determined by a proportion of infected tissue, then disease symptoms are developed, and the person is considered as infected and fall sick. If this critical level of virus concentration v c is such that v * 1 > v c > v * 2 , then people with weak immune response become ill while people with strong immune response remain healthy. The same property is preserved in the full model, where s 1 < s 2 leads to v * 1 > v * 2 for the same v 0 (see Appendix 2) .
However, if the initial virus load increases with the progression of epidemic, then v * 2 also increases and it can overpass the critical level. We consider this question in the next section.
We can now bring together the epidemiological model and the immune response model making a link between the number of infected individuals I, the initial viral load v 0 , and the number of susceptible individuals S. We suppose that the v 0 = v 0 (I) and S ≈ S(v 0 ) = S(v 0 (I)). Hence, the number of susceptible becomes a function of the number of infected. As before, let us assume that the total population consists of two sub-populations, one of them with weak immune response and another one with strong immune response. If the initial viral load v 0 is less than some critical valuev 0 , then the maximal viral concentration in the case of strong immune response remains less than the critical level v c which determines the appearance of the disease. This critical viral level can be determined from the equality v c =v 0 exp(W/s 2 ). Hence
where S 1 corresponds to the sub-population with weak immune response and S 2 to the total population. Under this assumption, the explicit dependence v 0 (I) is not needed. We suppose that the viral load is an increasing function such that v 0 (I) <v 0 for I < I c , and v 0 (I) ≥v 0 for I ≥ I c . Then instead of equation (4.1) we can write:
We can now find the solution:
Here
Hence, the number of infected individuals grows with exponent α 1 till some critical time t c , and with a larger exponent α 2 for t > t c .
Quarantine model. Another epidemiological model was introduced in [16] , in order to take into account removal of infected individuals to quarantine when the incubation period is finished and they manifest disease symptoms:
Here I(t) is the density of latently infected individuals (during the incubation period), τ is the length of the incubation period. In this case we have three time intervals: 0 < t < t c where t c is such that I(t c ) = I c ; t c < t < t c + τ , and t > t c + τ . We assume S(t) = S 1 , S(t − τ ) = S 1 for 0 < t < t c , and the solution of equation (4.3) is as follows:
where I 0 is the initial number of latently infected individuals, and µ 1 satisfies the equation
. This equation has a positive solution if the basic reproduction number R 1 = λS 1 τ is greater than 1. Let us recall that R 1 characterizes the number of newly infected individuals during the incubation period. The value t c can be found from the equality I(t c ) = I c , t c = 1 µ 1 ln Ic I 0 . During the second time interval t c < t < t c + τ , we have S(t) = S 2 and S(t − τ ) = S 1 . The solution of equation (4.3) writes:
where µ 2 satisfies the equality µ 2 = λS 2 (1 − θe −µ 2 τ ), θ = S 1 /S 2 . Since θ < 1, this equation has a positive solution and hence µ 2 > 0.
If t > t c + τ , then S(t) = S 2 , S(t − τ ) = S 2 ,
where I 1 = I c e µ 2 τ , and µ 3 satisfies the equation µ 3 = kS 2 (1 − e −µ 3 τ ). This equation has a positive solution if the basic reproduction number R 2 = λS 2 τ is larger than 1. Let us note that µ 1 < µ 3 < µ 2 . Therefore, the final growth rate is larger than the initial one but less than the intermediate.
Numerous data on seasonal influenza epidemic shows that growth of the number of infected individuals can occur in two stages. Both of them are characterized by an exponential growth rate but the second exponent is larger than the first one. Growth accelerates when the number of infected reaches certain threshold. It should be noted that usually several virus strains are observed during influenza epidemic. All of them show similar growth pattern.
Such dynamics can not be captured by conventional epidemiological models like SIR (susceptible-infected-recovered). These models predict exponential growth at the beginning of epidemic when the number of susceptible S can be considered as approximately constant or too large compared to initial number of infected. During the disease progression, S decreases, as a consequence the growth rate of disease progression also decreases. The number of infected passes through a maximum and then decays to 0.
The goal of this work is to develop a minimal mathematical model which can describe accelerated growth of disease spread. Our main hypothesis is that this acceleration can occur if we take into account the presence of different sub-populations, namely, with weak and strong immune response. In the beginning of epidemic, only (or mainly) people with weak immune response are concerned while at the second stage people with strong immune response also become susceptible. This transition takes place when the number of infected individuals becomes significantly large, and an uninfected individual has enhanced probability of interactions with infected individuals. This increases the viral load received during a short period of time and leads to the disease propagation in the sub-population with strong immune response.
In order to analyze how the initial viral load influences the disease initiation and progression, we consider an immuno-epidemiological model describing the concentration of viruses and immune cells in the organism. The strength of immune response in the model is determined by the rate of multiplication of immune cells. We show that for the same initial viral load, an individual with a weak immune system falls ill while an individual with a strong immune system does not develop disease symptoms. However, if the initial viral load increases, then the individuals with strong immune response also start developing the disease.
The combination of immunological and epidemiological models allows the description of the interplay between the population level and individual level [3, 4, 10, 14, 15] . We use here this approach to study a possible acceleration of epidemic growth rate.
We describe a generic situation applicable to different types of epidemics. During the influenza epidemic, the critical level of incidences when the growth accelerates is from 50 to 150 per week over 100 000 of the population (Figures 1, 3, 4) . Taking into account that the disease duration is about one week, then we get an estimate of the number of simultaneously infected individuals. The situation is different for the coronavirus infection because the disease duration is longer (about 4 weeks) and because the patients with manifested symptoms are put in quarantine. Therefore, we need to take into account here only the latently infected individuals during the incubation period. For a country with the total population 10 million and duration of the incubation period is one week, the critical level of epidemic is about 1400 new daily cases. This estimate does not take into account the restrictions on the population movement (it does not affect the influenza epidemic).
The immuno-epidemiological model developed in this work is not pathogen specific. We suppose that it can be also applied to coronavirus infection since the strength of adaptive immune response plays important role in the disease progression [8] . It should be noted that in some cases excessive immune response can damage lung tissue and provoke important clinical consequences.
Limitations of the model. The main hypothesis of this work is that epidemic progression begins with a sub-population characterized by weak immune response and, at a later stage, it continues in the whole population. Though this assumption seems plausible, we do not have direct confirmation that this heterogeneity of the population plays a significant role during epidemic spread.
A particular choice of the epidemiological model is not essential here since all of them give exponential growth for the number of infected individuals. The exponential growth rate is proportional to the density of susceptible individuals. Once we suppose that effective density of susceptible can increase during epidemic because different sub-populations are involved one after another, we can obtain the increase of the growth rate.
The immunological model used in this study is very simplified. It does not take into account for the involvement of many different cell types, intracellular regulation, time delay in virus production and clonal expansion of immune cells, and some other relevant aspects. However, it captures the main features of the interaction of viral infection with the immune response. Its simplicity allows to couple it with the epidemiological model.
Let us note that system (3.1), (3.2) can also be explicitly solved if the functions p(v) and h(v) are linear (and the second assumption is satisfied). The strength of immune response is characterized by the parameter p 1 which determines the rate of clonal expansion of immune cells in response to the antigen. A lesser value of p 1 Figure 3 : Data on seasonal influenza in UK [12] . The left graph shows the number of influenza-positive specimens from non-sentinel sources (blue line) and its approximation by an exponential (red line). The right graph shows the same curves in the logarithmic scale. Figure 4 : Data on seasonal influenza in Belgium [13] . The left graph shows the number of influenza-positive specimens from non-sentinel sources (blue line) and its approximation by an exponential (red line). The right graph shows the same curves in the logarithmic scale. corresponds to weaker immune response ( Figure 5 ), a larger value of p 1 to stronger immune response ( Figure 6 ). There are 7 trajectories in each figure starting from the same (pairwise) initial condition for some small values v = v 0 (initial viral load) and c = c 0 (initial concentration of antigen specific cells). If the initial viral load is sufficiently small, then in both cases the trajectories converge to the virus-free stationary point corresponding to the curing. However, in the case of strong immune response the maximal and the total virus concentrations are less. Therefore, strong immune response can prevent the development of the disease while weak immune response may not be sufficient for this. If the initial viral load is sufficiently large, then weak immune response can lead to a chronic infection or to death while strong immune response leads to virus eradication.