key: cord-0729237-307mrq0y authors: Berec, L.; Smycka, J.; Levinsky, R.; Hromadkova, E.; Soltes, M.; Slerka, J.; Tucek, V.; Trnka, J.; Smid, M.; Zajicek, M.; Diviak, T.; Neruda, R.; Vidnerova, P. title: The COVID-19 epidemic in the Czech Republic: retrospective analysis of measures (not) implemented during the spring first wave date: 2020-11-10 journal: medRxiv : the preprint server for health sciences DOI: 10.1101/2020.11.06.20227330 sha: 83967a33d9d9205d06094851622eed9d03d8eb81 doc_id: 729237 cord_uid: 307mrq0y The Czech Republic (or Czechia) is facing the second wave of COVID-19 epidemic, with the rate of growth in the number of confirmed cases (among) the highest in Europe. Learning from the spring first wave, when many countries implemented interventions that effectively stopped national economics (i.e., a form of lockdown), political representations are now unwilling to do that again, at least until really necessary. Therefore, it is necessary to look back and assess efficiency of each of the first wave restrictions, so that interventions can now be more finely tuned. We develop an age-structured model of COVID-19 epidemic, distinguish several types of contact, and divide the population into 206 counties. We calibrate the model by sociological and population movement data and use it to analyze the first wave of COVID-19 epidemic in Czechia, through assessing effects of applied restrictions as well as exploring functionality of alternative intervention schemes that were discussed later. To harness various sources of uncertainty in our input data, we apply the Approximate Bayesian Computation framework. We found that (1) personal protective measures as face masks and increased hygiene are more effective than reducing contacts, (2) delaying the lockdown by four days led to twice more confirmed cases, (3) implementing personal protection and effective testing as early as possible is a priority, and (4) tracing and quarantine or just local lockdowns can effectively compensate for any global lockdown if the numbers of confirmed cases not exceedingly high. could not be kept forever, as it does not only dive economics, but also negatively affects mental state of people. For example, an increase in divorce rate was reported in Wuhan, China, after quite stringent restrictions were eventually alleviated there (19) . Therefore, it is necessary to look back and assess effectiveness of the particular interventions adopted during the spring lockdown. Moreover, some other measures such as tracing and quarantining recent contacts of confirmed cases or implementing a geographically restricted lockdown were later applied or considered, and it is useful to look also at what might happen if they were alternatively implemented in the early phases of epidemic. Here we develop a detailed compartmental model to (i) reproduce the initial phase and lockdown period of COVID-19 epidemic in the Czech Republic, accounting for all implemented interventions, (ii) conduct a retrospective analysis of those interventions, and (iii) assess potential effects of alternative strategies proposed during the first wave or later. We structure our model by age, type of inter-individual contacts and space, thus allowing to explore and compare relative efficiency of a number of realistic intervention strategies that can be incorporated and timed explicitly (Methods). It contains a core epidemic layer, hospital layer, quarantine layer, and testing layer. In our model, we integrate public health data on the first wave in Czechia, results of sociological surveys performed before and during the first wave, real-time population movement data, and epidemiological parameters estimated from both public health and published data. To account for all these sources of information and their uncertainty, we use an Approximate Bayesian Computation framework based on massive super-computing simulations; this allows us to assess credibility intervals of our results given up-to-date knowledge on COVID-19 (Methods). To reach the goals we outline above, we first focus on the initial phase of epidemic and subsequent lockdown, spanning the period between 1 February 2020 and 7 May 2020. The first three cases of the COVID-19 epidemic in the Czech Republic were reported on March 1, 2020. Population-wide interventions soon followed (Table 1) . Till May 7, the restriction scheme was clear, there were no unreported local interventions, and the numbers of confirmed cases were high enough to match assumptions of our population model. Therefore, we use 1 February 2020 as our arbitrary starting date to cover a real epidemic beginning (the first three cases apparently got infected and showed symptoms before March 1 and some undetected asymptomatic cases were certainly imported, too) and 7 May 2020 as a date until which all those interventions apply. Based on the age-specific cumulative numbers of confirmed cases from this period, we infer uncertain epidemiological parameters (i.e., calibrate our model) and use the resulting setup as the baseline state for our retrospective analyses. The mechanistic character of our model allows for exact implementation of specific dates of intervention initiation (and relaxation), and for setting up factors that reflect impacts of these interventions on epidemiological as well as behavioral parameters A) B) Figure 1 : Sociological data used in our study. (A) Behavioral responses before and during the lockdown in Czechia, based on public opinion polls organized by the PAQ Reseach agency (Methods; www.paqresearch.cz) every week (in between the subsequent polls, data are linearly interpolated). The orange line is a proportional reduction in the number of contacts with respect to the pre-pandemic state (value 1), the red and blue lines are proportions of respondents that reported using face masks and increased hygiene. (B) Reduction of the weekly number of contacts in different environments due to the COVID-19 epidemic, according to the public opinion polls organized by the Median agency (Methods; www.median.eu). The respective proportional reductions in contacts that we use during the lockdown are thus 0.44 for home, 0.45 for work, 0.03 for schools, and 0.35 for other types of contact; our model distinguishes just these four types of contact (Methods). our model contains. Where possible, we base these factors on results of sociological surveys on behavioral responses before and during the lockdown (Methods, Fig. 1 ). We use the Approximate Bayesian Computation (ABC), a technique used to estimate parameters of complex models in genomics and other biological disciplines (1, 6) , including epidemic modeling (2, 17) , to integrate all our data sources, including their uncertainty. More specifically, using a high performance computer, we simulated 100,000 realizations of epidemic dynamics from the model, using parameter values randomly generated from prior distributions based on public health data in Czechia or literature review (Methods; blue distribution in Fig. 2A ). Subsequently, we used a rejection-sampling algorithm to select 0.1% of model realizations (100) that best matched the observed age-specific cumulative numbers of confirmed cases (Fig. 2B) . The parameters used to generate these selected realizations form an outer estimate of Bayesian posterior distribution of parameter values (1, 6) (red distribution in Fig. 2A ). The distributions of parameter values corresponding to selected simulation runs thus possess residual uncertainty in the parameters after the model is confronted with the observed data (1, 6) . We may think of the individual accepted parameter sets as different worlds A) B) Figure 2 : Approximate Bayesian Computation (ABC). (A) Prior parameter distributions based on available information are first built for every parameter to be estimated (blue distribution). (B) Parameter sets randomly generated from these priors are then used to run the model, and the sets for which the corresponding simulations provided a sufficiently good fit to the observed data were not rejected and were used to built posterior distributions for the estimated parameters; red distribution of panel A is the marginal distributions calculated from the resulting multidimensional posterior distribution. our observations of which are compatible with real observations, and we can study the epidemic of COVID-19 in any of these worlds or in all of them simultaneously (Methods). As Figs 2B and 3 indicate, our model is able to grasp dynamics beyond the COVID-19 epidemic in Czechia. We emphasize that while reasonably fitting the numbers of confirmed cases in any age cohort (Fig. 3 ) means reasonable fit of the total number of confirmed cases (Fig. 2B) , the converse is not true (not shown). Therefore, we think that age-structured models should always be calibrated on age-structured data, simply to be able assess effects of interventions that are often targeted to specific age cohort. A lot of epidemiological models developed to describe the COVID-19 epidemic have been prospective and were not calibrated on real observations, providing qualitative predictions forward in time (8, 10, 11) . Another large group of models are statistical in nature, attempting to retrospectively analyze an effect of applied interventions right from observed data (4, 9) . Our model in a sense bridges these two groups of models, providing a mechanistic description of the epidemic, enabling freely switch interventions on and off and modify their intensity, while at the same time being calibrated on a robust set of observed data. We consider the following scenarios that we are going to compare with the baseline one 5 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint Figure 3 : Baseline scenario. Temporal dynamics of the age-structured and total cumulative numbers of confirmed cases for the baseline scenario, covering the initial phase of epidemic and the following lockdown, for children (age cohort 0-19, A), adults (age cohort 20-64, B), seniors (age cohort 65+, C), and the probability density functions for the number of confirmed cases by 7 May 2020 (D). In panels A-C, blue curves represent real observations, while the orange lines are results of model simulations run for 100 selected posterior parameter sets. In panel D, the vertical lines represent the numbers of confirmed individuals by May 7, 2020. 6 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint (more detailed description of each scenario follows): R1 Lockdown interventions adopted four days before their actual date R2 Lockdown interventions adopted two days before their actual date R3 Lockdown interventions adopted two days after their actual date R4 Lockdown interventions adopted four days after their actual date R5 No contact reduction assumed, degree of personal protection as observed R6 Contact reduction as observed, no personal protection assumed R7 Dates of imposing protective measures and closing schools plus recommending home office swapped R8 Testing in March assumed as effective as later in April and May R9 Tracing and quarantine applied since April 1, 2020 R10 Local lockdowns applied since April 1, 2020 All these scenarios work with posterior parameter sets identified above; in each of the 100 selected posterior sets a change in policy is made and a corresponding model is run. An absolute impact of each scenario is examined, as well as its relative impact with respect to the baseline scenario (Methods). Experience from many countries demonstrates that one of the most important issues is timing of interventions; every single day of postponing them or of implementing them earlier counts (we are of course aware of other than epidemiological factors that enter the stakeholders' decision process, but we focus just on epidemiology in this study). Indeed, establishing every single intervention two or four days earlier or later produces significant differences in the number of confirmed cases by 7 May 2020. (Figs 4 and 5) . Thus, every day counts in the rate of responding to the incipient epidemic. As we show below, delay of four days in implementing lockdown roughly doubles the number of confirmed cases by May 7. The scenario R5 is inspired by a highly discussed issue of whether leave schools open or rather close them (4, 21, 14) . Since leaving schools open means in our model on average nearly 5 child-child contacts more per day, we extend the commonly assumed scenario of closing schools to asking more generally whether personal protection itself would suffice, without limiting any type of contact in the population. Not limiting contacts worsens the epidemic about five times (Fig. 6A ). On the other hand, keeping the amount of contacts limited as observed, but without any personal protection, produces twenty times more perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. sensitive to changes in contact structure than to changes in the factor reducing infection transmission due to protective measures (face masks, increased hygiene, interpersonal distance). The art of responding to any epidemic is not only in setting up the right restrictions and setting up them in time, but also setting them in an optimal order or implementing them as effectively as possible. With regard to our above finding that personal protection measures in the form of wearing face masks, using hand disinfection and keeping interpersonal distance appears more effective than just limiting social contacts, we suggest that implementing the former as early as possible is more effective (Fig. 7A) . Similarly, any decrease in the delay related to testing symptomatic cases weakens the epidemic spread (Fig. 7B ). Effects of two scenarios that have been used or considered to compensate an expected increase in the number of confirmed cases when lockdown was relaxed are now examined. We consider either tracing contacts of confirmed cases and their quarantining or local lock-10 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; down when a given threshold in the number of new cases per 100,000 inhabitants during the last seven days is exceeded in any specific county. Either intervention is implemented since the beginning of epidemic. Both scenarios suggest that both tracing and quarantine and only local lockdowns can successfully replace any global lockdown if intensive efficient enough (Fig. 8) ; no contact limitation was imposed in tracing and quarantine scenario and in counties not under local lockdown, while the degree of personal protection remained as under global lockdown. More quantitative way of assessing the effect of (not) implementing an intervention is to calculate effects of an intervention relative to the baseline (i.e., general lockdown) scenario when both simulations come from respective worlds (Methods; Fig. 9 . Obviously, intervening as soon as possible, sequencing the involved restrictions optimally and testing efficiently is by far a preferred strategy. We show that every four days of delaying the interventions double the number of confirmed cases by the end of period we follow, that personal protective measures are more efficient than reducing the number of contacts and that tracing and quarantine or just local lockdowns may compensate for relaxation, but only if some measures still remain in effect. These conclusions may seem trivial, and they actually are. The advantage of modeling epidemics is in providing some quantitative assessment of what intensity of interventions is not enough and which may hopefully suffice. We developed a mathematical model to analyze the spring first wave of COVID-19 in the Czech Republic. The model allows decoupling of protective as opposed to contact-limiting measures, setting up dates of implementing and relaxing of various types of interventions, and structures the population by age (three age cohorts: children, adults and seniors), type of contact (four types: at home, at school, in work and within community) and space (206 counties). Moreover, we included tracing and quarantining contacts of confirmed individuals and calibrated the model to the cumulative number of confirmed cases (so we needed a testing layer that modeled transition from real numbers of cases to the observed ones). We found that (1) personal protective measures as face masks and increased hygiene are more effective than reducing contacts, (2) delaying the lockdown by four days led to twice more confirmed cases, (3) implementing personal protection and effective testing as early as possible is a priority, and (4) tracing and quarantine or just local lockdowns can effectively compensate for any global lockdown if the numbers of confirmed cases not exceedingly high. Two types of models have generally been used to get an insight into the effects of nonpharmaceutical interventions in controlling COVID-19 epidemics. One of these groups, exemplified by Davies perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. Legend for the left panels as in Fig. 4 , the right panels show the proportions of counties (out of 206) that are closed at a given date. 13 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. 14 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint use parameter values found in the literature, often coming from different parts of the world, and are not calibrated on real observations in a given country. Another group of models, exemplified by Brauner et al. (4), Dehning et al. (9) , are statistical in nature, attempt to retrospectively analyze an effect of applied interventions right from the observed data, and often analyze many countries at once. Our model in a sense bridges these two groups of models. It provide a mechanistic description of the epidemic, enables freely switch interventions on and off and modify their intensity, while is at the same time calibrated on a robust set of observed data from the Czech Republic. Hence, it can be used to make prospective predictions, but can also look back and retrospectively analyze an effect of interventions applied in the past. One of its unique features is a use of sociological data to quantify a degree of compliance of established interventions, that is, a degree of contact limitation in various environments as well as a degree of compliance of obeying required protective measures (wearing face masks, using disinfection, keeping interpersonal distance). Another unique feature of our model is how it deals with uncertainty. Although a deterministic model, it is calibrated via a stochastic technique called Approximate Bayesian Computation (1, 6) applied to model parameters least known from either the literature or data we have for the Czech Republic epidemic. In this way, we generate a number of posterior parameter sets (or loosely said different hidden worlds) that all more or less result in observations we make on the numbers of confirmed cases in any of the three age cohorts, and ask how in any of these world the observations would look like if our responsive strategies change. Moreover, this technique deals with uncertainty in the model structure; the estimated parameters need not always be biology realistic (even if we try to reflect current knowledge in the respective priors), but if behavior of a strategy relative to the baseline scenario demonstrates consistent behavior across different worlds, we may have an increased confidence in the actual effect, since even if worlds are different they demonstrate similar relative impacts. Why do we apply our model just to Czechia and not to a wider set of countries? There are several reasons for that. First, the adopted detail of our model requires that: we do not have sufficient data for other countries (such as sociological data on behavior during epidemics, moving patterns of population). But also, other countries generally differ in a number of things, such as school system, community structure, structure of households, health system, testing strategy, etc. Nevertheless, there is no country-specific feature in the model structure, so once properly parameterized it can easily be used to any country affected by COVID-19. An issue commonly discussed in both country-specific studies (8) and studies spanning more countries (4), is whether to leave schools open or to close them. Closing schools turned to be among the first interventions in many countries (4), presumably because of governmental pandemic plans prepared tailored to influenza. This issue, causing quarrels already 100 years back during the Spanish influenza pandemic, still persists and studies 15 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; of COVID-19 can also be divided into those that claim relatively small effect (8, 14, 21) and those that rank closing schools among the most important interventions (4, and references therein). We do not follow this line of inquiry, mainly because it makes no sense to single out this specific intervention. Closing school requires parents of small pupils to remain home and thus affects home and work contacts, while large pupils tend to gather in shopping malls and other public places, affecting contacts within a broader community (this is one of the things that may differ among countries). To prevent just this shift of contacts elsewhere one needs to implement broader set of interventions beyond just closing schools, which is exactly what happened in Czechia. Shopping malls were closed shortly after in spring, and open wi-fi signal has been switched there now. Therefore, although schools can be closed in our model without affecting the other types of contact, we examined whether limiting contacts at large was more or less effective than personal protective measures, showing tat the latter is what should always be enforced first. Interestingly, Brauner et al. (4) rank wearing masks among interventions that decreased reproductive number relatively less than closing schools and businesses. Brauner et al. (4) claimed that a probable reason for this conclusion was that wearing face masks was among the last interventions implemented during lockdown in many countries, including Czechia (interestingly, many other countries enforced wearing face masks just following Czechia, seeing that it apparently works, which after some time admitted also WHO). We showed that implementing wearing face masks, together with other protective measures, among the earliest interventions would prevent many from containing infection during the spring first wave. On the contrary, the prospective study of Ferguson et al. (11) is more in line with our findings: although not fitted to any time series data, the authors find much less pronounced effect of closing schools than of implementing social distancing measures. Western civilization appears to be not much obedient in adherence to guidelines. An indirect indication of this is a much severe second wave in Europe compared to many South-Asian countries, even those with democratic governments. Any kind of lockdown is always a displeasing state that clashes with our perception of freedom. However, the lockdown is usually the last attempt to get an epidemic under control, and it turns out that is in not often necessary if we effectively adhere to less severe restrictions that usually come first. We suggest that a combination of timely application of protective measures and somewhat limited contacts, effective testing and contact tracing, possibly couples with local lockdowns, is apparently an optimal strategy. This may sound trivial. Mathematical models are needed to say how much is not enough, but also how much is probably unnecessarily much. Here we provide a description of our mathematical model, structured by age, type of contact and space. The model consists of four layers: a core epidemic layer, observation 16 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; layer, quarantine layer, and hospital layer. Although we do not report results regarding hospitals and deaths in this study, due to low number of affected people during the first epidemic wave, for completeness we give description of the hospital layer here. We first describe unstructured versions of all layers. Their extensions to age, various types of contact, and space then follow. Our core epidemic model is a variant of the classic SEIR model. Due to contacts with infectious individuals, susceptible individuals (S) may become exposed (E), that is, infected but not yet infectious (the process of infection transmission is described below). The exposed individuals then become asymptomatic for the whole course of infection (I n , with probability 1 − p S ) or presymptomatic for just a short period of time before becoming symptomatic (I a , with probability p S ). The I a individuals later become symptomatic, reducing their contacts with others by the factor r C . We assume that a proportion p T of individuals that become symptomatic decide to undergo testing for the presence of SARS-CoV-2 (I s ). They either contact their general practitioner and are sent for testing or go and pay by themselves. On the contrary, the proportion 1 − p T of symptomatic individuals (most likely those with very mild symptoms) decide not to undergo testing and rather stay at home (I h ). The I n , I s and I h individuals then recover (R). Since in the Czech Republic, deaths attributed to COVID-19 did not generally occur outside hospitals, we do not consider deaths in the core epidemic layer, but only in the hospital layer described below. Deaths outside hospitals were an important issue in countries heavily hit by the COVID-19 spring wave where hospital capacities were soon depleted. Considering discrete time, with one time step corresponding to one day, our core epidemic model consists of the system of six equations: The hitherto unexplained variable L[t] accounts for the imported COVID-19 cases from abroad, mostly from Italy and Austria. A list of all confirmed (symptomatic) imported cases is available at onemocneni-aktualne.mzcr.cz/covid-19. However, we do not introduce such imported cases as symptomatic. We assume they came earlier as exposed, and introduce them before they were actually tested positive (to account for delay between exposition and confirmation). Moreover, to account for the likely situation that some of the imported cases remained undetected as being asymptomatic for the whole course of infection, we divide the number of known imported cased by p S , the probability of exposed individuals eventually becoming symptomatic. The force of infection λ in the model (1) sums contributions from all infectious classes, that is, I n , I a , I s , and I h : Here, β is the probability of infection transmission upon contact between susceptible and infectious individuals, C is the contact rate (the mean number of other individuals an individual has an effective contact with per day), r β is a factor reducing the infection transmission probability for an asymptomatic individual relative to a symptomatic one, r C is a factor reducing the contact rate of a symptomatic individual relative to an asymptomatic one (having symptoms should force an individual to reduce contacts with others), and N [t] is the total population size at time t. The remaining model parameters, σ, ξ, γ s and γ n , represent probabilities at which individuals leave the respective model classes. They are related to the mean duration an individual spends in each such class (Table 7) , which in turn is defined and set in Table 6 . All model variables are summarized in Table 5 . We emphasize here that the calculation of λ will change (be expanded) further in the text, as we extend our model description. We assume that the period from the onset of symptoms, through sampling and subsequent processing, up to infection confirmation and case isolation is assumed to take d T days (testing delay). The length of this period was found to decrease during the course of epidemic as all involved steps were performed more efficiently. The infectious individuals are assumed always tested positive (we do not assume false negatives). We thus have equation for I s redefined as: where η = 1 − exp(−1/d T ) is the testing rate. The number of newly positively tested symptomatic individuals at day t thus equals η[t] I s [t]. Therefore, the total number of positively tested symptomatic individuals yet to be reported (B) is where κ is a (also possibly time-varying) publication rate. This rate is calculated as κ = 1− exp(−1/d P ), where d P is the period from case confirmation to case reporting (publication delay). The total number of reported positively tested symptomatic individuals (K) is therefore This simple testing procedure, assuming that all symptomatic individuals are at any time t tested and reported at the same rate, links the real number of symptomatic individuals to the number of reported positively tested symptomatic individuals, and requires two (possibly time-varying) parameters d T and d P ( A proportion p H of positively tested individuals (those with relatively severe symptoms) needs hospitalization. The remaining proportion 1 − p H of positively tested individuals (those that have only mild symptoms) are sent home to stay isolated (I z ) until recovery. Hence, we add a new equation A scheme showing structure of the model involving all elements defined up to now is sketched in Fig. 10 . Hospitalized individuals may follow several pathways, depending on the number of hospital states one considers. Whereas most published models have not consider a class of hospitalized individuals (15) , the others considered one to three hospital states: one to cover all hospitalized individuals (12) , two to distinguish between common hospital beds and ICUs (10) , and three to further detach ICU patients that need lung ventilators or ECMOs (23). We consider one hospital state, for which we need to introduce three parameters: the proportion p D of hospitalized individuals that eventually die, mean duration from hospital admission to death d HD , and mean duration from hospital admission to recovery d HR (Table 3) . We thus introduce two hospital classes, H D and H R , representing the hospitalized individuals that eventually die or recover, respectively, and compose the following 19 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint dynamic equations for these two hospital classes: Here α HX = 1 − exp(−1/d HX ) for X = D, R. Also, we need to introduce an equation for dead individuals and modify the equation for those that recover: The force of infection λ for the S individuals (2), now needs another modification: N [t] in its denominator is now the total population size at time t except those individuals that have already died (D), are isolated at home (I z ), or are hospitalized (H D and H R ). These four classes are expected not to transmit infection to others. Quarantine is ordered to anyone that has been traced as having a recent (epidemiologically relevant) contact with a positively tested individual. are quarantined on the same day an individual is positively tested for SARS-CoV-2). So, to find the number of contacts traced and quarantined at time t we need to start with individuals that were tested positively at time We assume that d B + 1 is the number of days we look back for the contacts of a positively tested individual. On any single day from t − d I to t − d I − d B , there is C actual contacts with others. We assume that a proportion Φ of the actual contacts are successfully traced (successful tracing ratio; Φ = 1 means all recent contacts are traced). Therefore, there is overall Due to quarantine delay, an epidemiological status of a contact of a positively tested individual may change between the moments of contact and quarantine. We assume that the successfully traced contacts are only from the classes S, E, I a , I n , I h and R, and are initially proportionally distributed among these, with the class I h however weighted by r C (mobility reduction due to presence of symptoms . We follow the epidemiological status of any such individual in time. In particular, any of the S individuals selected as contacts becomes infected during the contact with a positively tested individual with probability β when the latter is symptomatic and with probability r β β when it is asymptomatic. All traced individuals actually undergo a short separate epidemic that runs for d I + w time steps, where w goes from 0 to d B , with infection transmission applicable only in the first time step: , , , where we use the time index τ instead of t not to mix it with the main system of equations. All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; then go to the respective classes Q S , Q E , Q Ia , Q In , Q I h , respectively. Consideration of the quarantine layer requires modification of the force of infection. The eventual form thus is: (11) Moreover, the quarantined exposed individuals Q E become Q In when asymptomatic for the whole course of infection or Q Ia when eventually symptomatic. The Q Ia then either go to I s with probability p T or to Q I h with probability 1 − p T . The force of infection by quarantined individuals is reduced by a factor r CQ (due to a reduction in their contact rate). Denoting by α the rate at which individuals sent to quarantine as S (hence to class Q S ) return to the class S if not infected in the meantime (in which case they would go to Q E ), the equations related to quarantine or containing a quarantine class eventually are: As SARS-CoV-2 is known to differently impact children, adults and seniors (8), we distinguish three major age classes: 0-19 years (children), 20-64 years (adults), and 65+ years (seniors). These classes interact via the force of infection: infectious individuals of any age classes contribute to the force of infection of a susceptible individual of any age class. Both the probability of infection transmission upon contact β and daily number of contacts C are thus now 3 × 3 matrices, referred to below as the transmission matrix and contact matrix, respectively. To account for this differential impact, we assume that the elements of the transmission matrix are not identical yet rather β has the following structure: where β 1 is a transmission probability between two children, β 2 is a transmission probability between children and adults or seniors, β 3 is a transmission probability between two adults, β 4 is a transmission probability between seniors and adults, and β 5 is a transmission probability between seniors. We estimate parameters β 1 , β 2 , β 3 , β 4 , and β 5 by fitting our model to age-specific data on the cumulative number of confirmed cases (see below). In addition to modeling age-specificity of COVID-19, having β as a 3 × 3 matrix allows for considering interventions that may impact various of its elements differently. The contact matrix C describes the mean number of other individuals of any age cohort that an individual of an age cohort meets per day. In these matrices, columns represent those that find, while rows those that are contacted. We exploit this division when defining and exploring impacts of various realistic intervention strategies. Once infected, individuals of each age cohort proceed independently of individuals of other age cohorts (except quarantine, where we trace contacts to a positively tested individual that are from any age cohort). Also, we formally treat most model parameters as agespecific, with the same value for each age cohort if age-specific information is not available. Finally, since we know age of any imported COVID-19 case in the initial phase of epidemic, we assign each such case to the appropriate age cohort. into well-defined 206 counties. For each county, the population size and its distribution into the above-defined three age classes is known (Czech Statistical Office, www.czso.cz/ csu/czso/home). Moreover, we compose and use a 206 × 206 mobility matrix that gives mean daily mobility patterns of individuals between all pairs of those counties (proportions of individuals travelling per day from a county to another). Only individuals from the classes S, E, I a , I n , and R are allowed to travel in our model, unless interventions are implemented to restrict their movement, too (see the section on interventions below). At each time step (one day), our age-structured model is first run (independently) in each county, and then the mobility matrix is applied to the updated county populations. Regarding the imported cases in the initial phase of epidemic, we have information about the broader region in the Czech Republic each case comes from, so we assign a random county from the corresponding region to each such case. The mobility matrix was constructed by averaging mobility patterns obtained from a telecommunication company across two weeks. Two such matrices were used, one representing the normal state using data from January 2020, and the other representing the lockdown state using data from the second half of March 2020. Both matrices were recalculated to account for the company's market share. Furthermore, movement matrices were adjusted to account for average visiting time of 6 hours, i.e. all movement intensities were divided by 4. Due to the lack of age-specific data, these matrices were assumed identical for all three age classes. Several scenarios for intervention implementation and relaxation are tested for their effect on the course of epidemic in the Czech Republic. Our baseline scenario considers all interventions that were in effect during lockdown established during March 2020 (Table 1 ). In modeling the various interventions, we exploit a division of the contact matrix C into four matrices describing contacts at home (C H ), school (C S ), work (C W ), and other types of contacts (C O ). Starting from the appropriate dates listed in Table 1 , we multiply the respective matrices by factors 0.44 for home, 0.45 for work, 0.03 for schools, and 0.35 for other types of contact. Moreover, personal protection, activated on 19 March 2020, including wide use of disinfection, wearing face masks on public, and keeping interindividual distances of more than 2 metres on public, was modeled so that all elements of the transmission matrix β were reduced by a given factor, corresponding to 88% efficiency (averaging masks and hygiene and computing mean over the high efficiency data). All these numbers are based on results of public opinion polls organized by two agencies during lockdown: the PAQ Reseach agency (www.paqresearch.cz) and the Median agency (www. median.eu). Results of these polls, summarized to our modeling purposes, are provided in Fig. 1 . These data show that during the second half of March and essentially the whole April, contacts of all kinds have been largely reduced while personal protection increased. All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint Values of several model parameters will always remain uncertain, of which the transmission matrix β is commonly the principal of those. This and some other model parameters, listed in Table 4 , are estimated by fitting the variable K (representing the cumulative number of confirmed cases in any age class and county), summed over counties but not age, to the age-specific time series on the reported cumulative numbers of confirmed cases in the Czech Republic. There are many ways to meaningfully perform model calibration on real-world data and many optimization and filtering methods exist (24). We use the Approximate Bayesian Computation (ABC), a technique used to estimate parameters of complex Markov models in genomics and other biological disciplines (1, 6), including epidemic modeling (2, 17) . The major advantage of this method is that it naturally works with all sources of uncertainty acknowledged in the model. At the same time, ABC does not rely of likelihood calculations and in case of sufficient computation power it can be used with models of virtually any complexity. The ABC procedure with rejection sampling we used here consisted of three steps. First, we used our model to simulate summary statistics for calibration (number of detected cases in age categories, state variable K) N times (100,000), drawing the uncertain model parameters from prior distributions based on literature (see above); selected prior distributions for estimated parameters are also given in Table 4 . Second, we compared the simulated summary statistics with the observed one, using the Euclidean distance D. Third, we selected model simulations that satisfy D < , where was chosen to pass 0.05% of simulations into the selected set. Given the used summary statistics are informative, the distribution of parameter used for such selected simulations is known to converge from outside to the Bayesian posterior distribution of parameter values with going to 0, and is referred to as approximate posterior (1) . The choice of and N in ABC is driven by compromise between computation power, smoothness and accuracy of approximate posterior. The set of accepted sets of parameters thus allows us to evaluate remaining parameter uncertainty, given the available data and adopted summary statistics (1, 6) . This is crucial to realize, since although different parameter sets may similarly fit the available data (have similar summary statistics), and often provide similar short-term predictions, they may demonstrate significant differences in longer-term predictions and in the interplay with intervention policies. To accommodate this uncertainty, we do not evaluate only an absolute impact of implementation or relaxation of an intervention. We also present its impact relative to the baseline (here the actually adopted full lockdown) scenario. If such a relative impact is consistent over whole posterior distribution of parameters, we may have confidence in its potential effect. To apply the ABC technique, we use the abc package in R (7), modified to work with non-normalized summary statistics. 25 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint The prior distribution for the latent period (d E ) was set as normal with mean 5.08 days and standard deviation 0.5 days, based on He et al. (13) . The duration of the infectious period before the onset of symptoms was reported in the literature as ranging from 1-3 days (22) . The mean duration of test positivity for viral RNA (a proxy for infectiousness) was reported in the literature as 25.2 days for symptomatic patients (d S ) and and 22.6 days for asymptomatic patients (d N ) (18) . The prior probability of factor reducing the infection transmission in a/presymptomatic individuals (r β ) was set as uniform distribution from 0.4 to 1, respecting a number of studies. In particular, He et al. (13) re-analyzed results previously published in local Chinese journals to show a somewhat higher probability of virus transmission in symptomatic individuals (the asymptomatic infection rate to be 46% of the symptomatic one, with the 95% CI (18.48%-73.60%)), but the source of primary data could not be verified. The value of 0.55 was assumed by Domenico et al. (10) in setting up their model, based on Li et al. (16) but for "undocumented" rather than "asymptomatic" patients. On the other hand, Chen et al. (5) showed no statistically significant difference in the transmissibility of the virus between symptomatic and a/presymptomatic cases. The percentage of asymptomatic individuals out of those positively tested for COVID-19 (p S ) has been estimated by a number of studies and the estimates vary significantly as does the quality of the evidence. The range of the percentage of the asymptomatic individuals among those positively tested reported is 4.4-89%. As a result, we use the uniform distribution on the whole interval 0-100% as a prior. All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted November 10, 2020. ; https://doi.org/10.1101/2020.11.06.20227330 doi: medRxiv preprint IAG BSR06-11, OGHA 04-064, HHSN271201300071C) and from various national funding sources is gratefully acknowledged (see www.share-project. org). This work was supported by the Ministry of Education Approximate Bayesian Computation in evolution and ecology HIV with contact tracing: a case study in approximate Bayesian computation Data resource profile: the survey of health, ageing and retirement in Europe (SHARE) The effectiveness of eight non-pharmaceutical interventions against COVID-19 in 41 countries. medRxiv The epidemiological characteristics of infection in close contacts of COVID]-19 in Ningbo city Approximate Bayesian Computation (ABC) in practice abc: Tools for Approximate Bayesian Computation (ABC) Age-dependent effects in the transmission and control of COVID-19 epidemics Inferring change points in the spread of COVID-19 reveals the effectiveness of interventions Impact of lockdown on COVID-19 epidemic inÎle-de-France and possible exit strategies Impact of non-pharmaceutical interventions (NPIs) to reduce COVID-19 mortality and healthcare demand Modelling the COVID-19 epidemic and implementation of populationwide interventions in Italy Estimation of the basic reproduction number, average incubation time, asymptomatic infection rate, and case fatality rate for COVID-19: Meta-analysis and sensitivity analysis Was school closure effective in mitigating coronavirus disease 2019 (COVID-19)? time series analysis using Bayesian inference Projecting the transmission dynamics of SARS-CoV-2 through the postpandemic period Substantial undocumented infection facilitates the rapid dissemination of novel coronavirus (sarscov-2) The epidemiological fitness cost of drug resistance in Mycobacterium tuberculosis Asymptomatic infection and atypical manifestations of COVID-19: Comparison of viral shedding duration China's divorce spike is a warning to rest of locked-down world Projecting social contact matrices in 152 countries using contact surveys and demographic data School closure and management practices during coronavirus outbreaks including COVID-19: a rapid systematic review Presymptomatic transmission of SARS-CoV-2 -Singapore We thank The Institute of Health Information and Statistics of the Czech Republic (IHIS CR; www.uzis.cz) for providing us detailed data on many characteristics of confirmed cases of COVID-19. We also acknowledge the PAQ Reseach agency and in particular Daniel Prokop (www.paqresearch.cz) and the Median agency (www.median.eu) for sharing with us results from their public opinion polls. The mobility data were provided by the project "RODOS Transport Systems Development Centre" which is supported by the Technology Agency of the Czech Republic. VT was supported by the QuantiXLie Center of Excellence, a project cofinanced by the Croatian Government and European Union through the