Model Description¶
Basic model¶
We model time \(t\) as an integer with the meaning of “days since day 0”, where day 0 can be chosen at will.
The primary dynamic variable \(\delta N_{\mathrm{inf}}(t)\) is the number of people being infected on day \(t\). The symbol \(\delta\) indicates that this is an increment per day, not the total number of cases. The total number of people that were infected between times \(t_0\) and \(t_1\) is:
For a population without immunity, the mean number of people infected by a single person is, by definition, given by the symbol \(R_0\) (the basic reproduction number). For a non-zero number of immune people \(N_{\mathrm{imm}}\), the actual mean number of people infected by a single person \(R_e\) equals
where \(N_{\mathrm{pop}}\) is the total number of people in the population.
The mean time between successive infections is called the generation time \(t_g\). It means that, on average \(t_g\) days after being infected, a person infects \(R_e\) further individuals with the virus. In practice this “passing on the virus” is spread over a time period, with different probabilities for different times \(t\) after the day of infection \(t_i\). We describe this “time spreading function of reproduction” as \(P_r(t-t_i)\), with obeys the following normalization:
because \(P_r(t-t_i)\) is a probability density function (PDF). The mean generation time is then its first moment:
With \(P_r(t-t_i)\) and \(R_e(t)\) we can now construct the infectivity function \(r(t_i,t)\)
which is the mean number of persons being infected on day \(t\) by an individual who was infected on day \(t_i<t\). The \(t-1\) in the \(R_e(t-1)\) is because \(R_e\) is not yet computed for time \(t\) when it is used in the computation of \(r(t_i,t)\).
Note that the use of \(R_e(t-1)\) for the computation of the infectivity function is the only place where a non-linear feedback is involved.
In the basic version of our model we make the following simplified assumption: A person infected on day \(t_i\) is assumed to be able to infect other people during a period of days from \(g_0\) up to (and including) \(g_1\) days after day \(t_i\). During that period, in this simple model, the infectivity \(r\) is constant, while before and after, the infectivity \(r\) is zero:
and the generation time is simply \(t_g = (g_0+g_1)/2\).
With \(r(t_i,t)\) defined, the epidemic dynamic equation for \(\delta N_{\mathrm{inf}}(t)\) becomes
The system of equations can be closed by assuming, for instance, that all (and only those) people that have been infected are immune:
and by assuming \(N_{\mathrm{pop}}\) to be constant and given, e.g. \(N_{\mathrm{pop}}=83,000,000\) for Germany.
The integration of this epidemic dynamic equation is straightforward: time-stepping in steps of one day. For \(R_0>1\) this leads to an exponential growth of \(\delta N_{\mathrm{inf}}(t)\) that saturates as a result of herd immunity.
The dynamics of the epidemic is entirely described by \(\delta N_{\mathrm{inf}}(t)\). Other quantities, such as the current number of ill people, the number of people in the hospital and on an intensive care unit (ICU), the number of people that have recovered or died, can all be computed from \(\delta N_{\mathrm{inf}}(t)\).
Derived quantities¶
For practical purposes, an epidemic dynamic model should be able to predict how many people will, at some time \(t\) in the future, be in a particular “state” (e.g. having symptoms, being in the hospital, occupying an ICU bed, or having deceased). Averaged over many different progressions of the disease, one can define a probability function \(p_{\mathrm{state}}(t-t_i)\) for a person that has become infected on day \(t_i\) to be in that state on day \(t\). These probability functions (PF) are, however, not normalized, because they are not probability density functions (PDFs), and because not every individual will ever acquire such a state. However, by their nature as PFs, they must obey
at all \(t\ge t_i\), and \(p_{\mathrm{state}}(t-t_i)=0\) for \(t<t_i\). Typically the numbers one gets from databases or reports is the fraction of people who acquire a certain state, e.g. that about 1% of reported cases get admitted to an ICU. That information is not sufficient to estimate the number of required ICU beds. The mean duration \(T_{\mathrm{ICU}}\) of the stay at the ICU, the mean time after infection \(t_{\mathrm{ICU}}\) that person is admitted to an ICU, and the probability \(p_{\mathrm{ICU}}\) that the person will enter such a state are required to make the simplest possible model of \(p_{\mathrm{ICU}}(t-t_i)\). In general this simplest model would look as follows:
where “state” can be replaced by “ICU”, “ill”, “hospital”, “recovered” or any other state of interest. More realistic expressions for \(p_{\mathrm{state}}(t-t_i)\) can be generated, for instance, by convolving this function with a Gaussian to account for the variance in \(t_{\mathrm{ICU}}\) and \(T_{\mathrm{ICU}}\).
The computation of the number of people in that state at time \(t\) is then straightfoward:
Adding features to the basic model¶
It is relatively straightforward to add features to the basic model. Here we give a few examples.
Effect of sudden policy changes¶
The implementation of a (partial) lockdown on day \(t_{\mathrm{lock}}\) will suddenly change the \(R_0\). In EpiDemo this can easily be implemented. One asks the model to integrate up to \(t_{\mathrm{lock}}\). Then one changes the \(R_0\) variable, and asks the model to continue integrating until the end.
Effect of a gradually implemented vaccination campaign¶
Already built into EpiDemo is the possibility to prescribe an array of values for \(N_{\mathrm{vacc}}(t)\) for the entire modeling time domain, which represents the number of people that have been vaccinated up to that point. EpiDemo then calculates \(N_{\mathrm{imm}}(t)\) as follows:
Effect of waning immunity¶
There are several ways by which a waning natural immunity can be implemented. The simplest being an exponential decay: Instead of writing \(N_{\mathrm{imm}}(t) = N_{\mathrm{inf}}(-\infty,t)\) one can write:
with \(\tau\) being the immunity waning time scale.
Multi-group version of the model¶
A real society is heterogeneous. Certain communities stay mostly among each other, having little interaction with other communities. An infection within such a community will, in the ideal (and unrealistic) case of zero interaction with the outside world, develop its own dynamics, with its own \(R_0\) value. Since the dynamics is exponential, even a small difference with the value of \(R_0\) outside that group will, after some time, lead to vast differences in the fraction of infected people compared to the outside world.
Suppose a group is so isolated that it managed to avoid infection altogether, even though already 5% of the population outside that group has become infected. Suppose that the behavior in this group is such that contacts are much more close and common than in the population at large (perhaps due to a false sense of security). If, for whatever reason, a single individual within this group does get infected, the behavior of this group may set their value of \(R_0\) to be much larger than that of the outside world. As a result, the virus spreads like wildfire inside this group, and soon the fraction of infected within the group greatly exceeds that of the more careful population outside.
If that group is large enough, the total number of registered cases in the country may suddenly see a sharp rise seemingly coming out of nowhere, as the total reported new infections of the group briefly outperforms the rest of the country, before the group reaches herd immunity or changes its behavior.
Since a group is never perfectly isolated, once the number of cases in the group massively dominates that of the outside world, even a small level of cross infection to the outside can then dominate the outside infection rates as well.
In Germany the massive outbreak occurred at a slaughterhouse in the region of Gütersloh in the summer of 2020 is likely an example of this, where the virus could, over an extended period of time, multiply in a relatively isolated community until a large fraction of that community was infected (see e.g. here). In the RKI data it can be clearly seen as a “bump” in the total number of reported new infections around June, which then spread to other regions of Germany and dissipated as a result of the overall negative \(R_0-1\) in Germany at that time.
This shows that even if one group is relatively isolated from the rest of society, if an epidemic is out of control in this group, the epidemic dynamics of society can become “enslaved” by the epidemic dynamics within the semi-isolated group (see e.g. the Keeling and Rohani book).
Note that this is different from a “superspreading event”, in which a single (or a few) people infect a very large number of people all at once due to circumstances that are exceptionally conducive to viral spread.
Most “groups” within a society are not so strongly isolated, and if contacts between groups are very common, there is no need to treat the groups separately in a model. The mean \(R_0\) value will then be enough to make reliable predictions.
But there are intermediate cases that are of particular interest. For instance, the attempts to isolate the elderly. Particularly the elderly in nursing homes are, out of necessity, in regular contact with nursing personnel. The degree of care taken to avoid cross infection from nursing personnel to patients determines the rate of cross-infection. How careful do they have to be to be successful in case the epidemic gets out of hand? Another example is the role of subcommunities of young people, some of whom may behave less careful than the rest of the population. How dangerous is this? None of these questions can be easily quantitatively answered, because the knowledge of the parameters is unavailable. But some simple scenarios can be tested.
In the multi-group model, \(\delta N_{\mathrm{inf}}(t)\) is replaced by \(\delta N_{\mathrm{inf},k}(t)\), where \(k\) is the index of the group. Python counts array indices from 0, so \(k\;\in [0,n-1]\) with \(n\) being the number of groups. The epidemic dynamic equation then becomes
where \(r_{kl}(t_i,t)\) is the infectivity function for infection from group \(l\) to group \(k\):
where we assume that \(P_r(t-t_i)\) is the same for all groups. Note the order of the indices: \(k\leftarrow l\), which is chosen to conform to the standard matrix multiplication index order. The “reproduction matrix” can be computed from the “basic reproduction matrix” \(R_{0,kl}\)
The values of the matrix \(R_{0,kl}\) have to be specified. The sums of the columns
give the the basic reproduction number for the average individual within group \(l\), where \(R_{0,ll}\) is the infection within the group and \(R_{0,l}-R_{0,ll}\) the infection to the other groups.
The off-diagonal elements of \(R_{0,kl}\) are tricky, because they depend not only on the behavior of people within the group, but also on how many people from the other groups they interact with. That depends also on the sizes of the groups. The matrix is thus also not necessarily symmetric. If an individual from group \(k\) meets with an individual from group \(l\), and both equally much expose the other, then \(R_{0,kl}N_{\mathrm{pop},l} =R_{0,lk}N_{\mathrm{pop},k}\).
For the rest, the model functions exactly the same as the one-group model.