key: cord-0658012-1sr8gs62
authors: Kantner, Markus; Koprucki, Thomas
title: Beyond just"flattening the curve": Optimal control of epidemics with purely non-pharmaceutical interventions
date: 2020-04-20
journal: nan
DOI: nan
sha: 96f466bfcc2b326f6cb048261dc8cee08731b96e
doc_id: 658012
cord_uid: 1sr8gs62

When effective medical treatment and vaccination are not available, non-pharmaceutical interventions such as social distancing, home quarantine and far-reaching shutdown of public life are the only available strategies to prevent the spread of epidemics. Based on an extended SEIR model and continuous-time optimal control theory, in this paper the optimal non-pharmaceutical intervention strategy is presented for the case that a vaccine is never found and complete containment (eradication of the epidemic) is impossible. In this case, the optimal control must meet competing requirements: First, the minimization of disease-related deaths, and, second, the establishment of a sufficient degree of natural immunity at the end of the measures, in order to exclude a second wave. Moreover, the socio-economic costs of the intervention shall be kept at a minimum. The numerically computed optimal control strategy is a single-intervention scenario that goes beyond heuristically motivated interventions and simple"flattening of the curve". Careful analysis of the resulting evolution of the time-dependent effective reproduction number reveals, however, that the obtained solution is in fact a tightrope walk close to the stability boundary of the system, where socio-economic costs and the risk of a new outbreak must be constantly balanced against one another. The model system is calibrated to reproduce the initial exponential growth phase of the COVID-19 epidemic in Germany.

shutdown of entire territories, are the only available measures. The interventions are aimed at "flattening the curve", i.e., reducing the transmission rate of the infectious disease in order to break the exponential growth of the epidemic.

In the case of the currently spreading COVID-19 epidemic caused by the new SARS-CoV2 coronavirus [35, 36] , the fundamental concern of the mitigation measures is not to exceed the available number of intensive care unit (ICU) beds, in particular for respiratory support or extracorporeal membrane oxygenation, in order to prevent actually avoidable deaths [22] . Since the outbreak of the epidemic, a large number of simulation studies has been published that use mathematical models to simulate the progression of the disease in order to assess the potential of various intervention scenarios and to estimate the corresponding demands on the health care system [1, 2, 4-6, 8, 9, 12, 15, 19] . Moreover, mathematical models are employed to deduce important epidemiological parameters via advanced data assimilation methods [4, 7, 14, 28] .

As long as no vaccine is available, a possible intervention strategy could be aimed at keeping the number of simultaneously infected below a certain threshold [23, 33] . This threshold can be either a manageable upper limit for tracing programs or, in extreme cases, the maximum number of intensive care beds to prevent the collapse of the health care system. In order to achieve this goal, non-pharmaceutical measures must be precisely coordinated, both in terms of time and intensity, to prevent a sudden increase in the number of cases. Due to the inherently nonlinear dynamics of epidemics, the optimal control of infection transmission is a non-trivial problem and heuristic intervention strategies such as simple bang-bang controls [21, 30, 31] or cycling strategies are typically far from optimal [5, 34] . Moreover, this type of interventions also might include ineffective measures that should be avoided for socio-economic reasons.

The central scope of this paper is to compute the optimal transmission rate, i.e. the optimal evolution of the effective reproduction number, on the basis of an extended SEIR (susceptible-exposed-infectious-recovered) model and continuous-time optimal control theory following a variational principle [17] . The scenario considered here is the (hopefully unlikely) case in which an effective vaccine is impossible or never found and the epidemic must be controlled with purely non-pharmaceutical measures. Furthermore, the scenario excludes the possibility of complete containment ("eradication of the virus"), as this is beyond the applied modeling framework (macroscopic mean-field models). Then, optimal control must pursue competing objectives: On the one hand, the number of disease-related deaths must be minimized by strictly avoiding an overload of intensive care treatment capacities. On the other hand, however, in the long run, at the end of the measures, sufficient natural immunity must have been built up in the population to exclude a second outbreak of the epidemic. Furthermore, the socio-economic costs of the intervention shall be kept at a minimum. The corresponding optimal control obtained in this paper, is a single-intervention strategy, which keeps the number of simultaneously infected almost constant over a sufficiently long period of time by means of a precisely tuned temporal course of the effective reproduction number. The underlying mathematical model for the progression of the epidemic and the estimation of the demand for intensive care resources is described in Sec. 2. The optimal control problem is derived in Sec. 3 and the results are described in Sec. 4. One of the major findings of this paper is that the numerically computed optimal control obeys two fundamental stability criteria, which impose an upper limit on the effective reproduction number and its rate of change on the way out of an initial lockdown with minimal socio-economic costs. The parameter adjustment to reproduce the exponential growth phase of the COVID-19 epidemic in Germany is described in the Appendix. (1). The function u (t) describes a modification of the transmission dynamics due to non-pharmaceutical interventions.

(b) State-dependent mortality rate f as a function of the number of patients in a critical state requiring intensive care. The mortality rate grows rapidly if the number of critical patients exceeds the number of available ICUs C 0 . Inset: The solid line is the regularized mortality rate (2b) that is used in the computations throughout the paper.

Mathematical modeling of the spread of epidemics is an indispensable tool to project the outcome of an epidemic, estimate important epidemiological parameters by adjustment to observed data and to make predictions for different intervention scenarios. Compartment models [3, 11, 13] , where the population is divided into different macroscopic sub-populations, such as susceptible, infectious, recovered etc., are a simple but effective tool to model the progression of epidemics. In contrast to complex (but more realistic) stochastic agent-based models, deterministic mean-field models are limited to the description of the average infection dynamics in macroscopic (sub-)populations, but allow for fast parameter scans and a straightforward application of continuous-time optimal control theory [17] . In this paper, an extended SEIR model, similar to that proposed by Neher et al. [18] , is used to model the spread of an epidemic and to estimate the number of patients in a critical state that require intensive care. Similar models are described in Refs. [14, 34] . For the sake of simplicity, vital dynamics (except for disease-related deaths), seasonality effects and different age groups are neglected. The total population can be subdivided into susceptible S, exposed E, infectious I, hospitalized H (severely ill), critical C, recovered R (i.e., immune) and deceased D individuals. The model equations reaḋ

The group of susceptible persons (S) is vulnerable to infection through contact with infectious persons (I) who may transmit the disease to susceptible persons. The infection probability is determined by the transmission rate β and the number of infectious Fig. 2 (a) Evolution of the epidemic without interventions (u (t) ≡ 1) for initially 20 exposed individuals. The number of available ICUs was set to C 0 = 30 000. The inset shows the overflow in ICU demand, which leads during a period of about 57 days to an increased mortality rate according to Eq. (2). (b) Same as in (a) but on a logarithmic scale. The markers indicate the estimated number of cumulative cases (see Appendix C) and the reported numbers for ICU demand and deaths during the early phase of the COVID-19 epidemic in Germany. The first disease-related fatalities were reported on March 9, 2020 (day number 20 in the simulation). Social distancing measures, which came into force nationwide in mid-March, have flattened the initial exponential growth.

individuals. The newly infected persons (E) become infectious themselves only after a latency period γ −1 l (which must not be confused with the incubation time). Infectious individuals either recover or turn severely ill after an average period γ −1 i . Severely ill patients (H) can either deteriorate into a critical state (C) or recover after a period γ −1 h . The group of recovered persons (R) is assumed to be immune against new infections. Patients in a critical state either stabilize to the severely ill state or die from the disease on a time scale γ −1 c . Moreover, m is the share of infectious individuals that are asymptomatic or have at most mild symptoms, c is the fraction of severely ill patients that become critical and f is the fraction of critical patients that are going to die from the disease. Finally, the time-dependent function u (t) describes a modification of the transmission rate due to non-pharmaceutical interventions such as social distancing, hygiene measures, cancellation of mass events etc. Here, u = 1 means no intervention, and u = 0 corresponds to the extreme case of total isolation of the whole population.

The model system is illustrated in Fig. 1 (a) .

The disease-related mortality grows tremendously as soon as the number of critical patients exceeds the number C 0 of available ICUs. This is modeled by a state-dependent average fatality rate

As long as every critical patient can be served with an ICU (C ≤ C 0 ), the fatality rate is a constant f = f 0 . As soon as the ICU resources are exceeded, an increasing fraction of the critical patients dies with a higher rate f 1 > f 0 , which on average results in the state-dependent fatality rate (2a). Here, f 1 = 2f 0 is assumed. In the following, the 

is used, in order to avoid problems due to the nondifferentiability at C = C 0 . The function f (C/C 0 ) is shown in Fig. 1 (b) . The statedependent fatality rate proposed here is considered to be more realistic than the discontinuous model employed in Ref. [5] .

The basic reproduction number

can be thought of as the expected number of cases (without intervention) that is directly generated by one case in a population where all individuals are susceptible to infection. The effective reproduction number

depends on time and includes the impact of intervention measures. 

In the scenario outlined in Sec. 1, where an effective vaccine is never found, the optimal transmission control u (t) due to non-pharmaceutical interventions is required to avoid (i) ICU overflow (more patients in a critical state than available ICUs) but at the same time (ii) exclude a second wave of the epidemic after the end of the measures. The optimal solution is computed by minimizing the index functional

is the terminal cost function. The first term in Eq. (5b) describes the number of diseaserelated deaths D (T ) at the end of the epidemic, which should be minimized. As the increment of the disease-related deaths depends on the state-dependent fatality rate, see Eq. (1g), this condition implies that the ICU capacities must not be exceeded. The second term of Eq. (5b) controls the number of susceptible individuals S (T ) at the end of the epidemic. In order to approach a stable, disease-free stationary state, the share of susceptible individuals on the total population must be less than R −1 0 at the end of the intervention, see Appendix A. Here, the desired final state is chosen as S f = 0.95 N/R 0 , see Fig. 3 (b). The function

is convex on the whole domain x ∈ [0, ∞). It appears also in the last term of Eq. (5a)

as an intermediate cost function, which provides an abstract measure for the total socio-economic costs caused by the intervention. The term is minimal and zero if no intervention is applied C (1) = C (1) = 0, see Fig. 3 . The advantage of using (6) over the commonly used quadratic cost functions is that "unphysical" negative values of u are a priori excluded. The control parameter P balances between the competing objectives of minimal disease-related deaths (first term), while attaining at the same time a minimum number of cases to enforce an upper limit on S (T ) (second term). Simultaneously, the total (socio-economic) costs of the intervention (third term) shall be minimized. Ramping up P puts an increasing emphasis on minimizing the disease-related deaths.

The final time T must be taken sufficiently large (several years, depending on C 0 ) in order to avoid finite size effects. Since the objective is to steer the outcome of the epidemic into a stable, disease-free, stationary state in which the epidemic is receding, the results are practically independent of T (provided it is large enough).

Following Pontryagin's maximum principle [17] , the optimality condition reads

Finally, the co-state equations and the final time conditions are obtained aṡ

Together with the initial conditions x (0), the system (1), (9)-(10) represents a nonlinear two-point boundary value problem. The full set of equations is given in Appendix D. Numerical solutions are obtained by using Matlab's built-in routine bvp4c [29] in combination with an analytic Jacobian matrix and a step-size adaptive homotopy method, where the control parameter P is gradually ramped up while always using the result of the previous step as initialization. The procedure is initiated from the numerical solution of the initial value problem (1) without interventions, see Fig. 2 .

With optimal control of the transmission rate (in the sense of Sec. 3) via accordingly steered non-pharmaceutical interventions, the epidemic develops dramatically different from the uncontrolled case, see Fig. 4 .

In the optimized scenario, the measures described by u opt (t) begin with a strict "lockdown" that is built up over a period of about 25 days (starting around day 25), see Fig. 4 (a, b) . The effective reproduction number (4) must be held below one R eff < 1 for about 11 days. This strict initial intervention breaks the early exponential growth and damps the peak number of infected such that an overshoot of the patients in a critical condition is just barely avoided, see Fig. 4 (b). The shape of the initial lockdown therefore reflects the long-term goal of the scenario under study, namely the establishment of a just sufficient degree of natural immunity at the end of the measures. This objective requires that, depending on R 0 , an enormous minimum number of infections are passed through and treated over a long period of time without exceeding the capacity limits of the health care system ("flattening the curve"). Remember that the possibility of vaccines is excluded in this scenario. The initial lockdown is followed by a longer period (about 300 days), during which the intensive care system is constantly stressed by slightly less than C 0 patients in a critical condition (here C 0 = 30 000 was assumed), 

In the optimized scenario, the number of susceptible individuals S (T ) /N is reduced slightly below the critical value R −1 0 , which guarantees herd immunity and rules out a second wave of the epidemic. (b) Zoom on the lower part of (a). The optimal transmission control ensures that the available number of ICUs is not exceeded by patients in a critical state: C (t) < C 0 for all t ∈ [0, ∞). A more detailed plot of the ICU load is given in Fig. 5 (c) . in the deterministic model (1) at all. During this period, which will be denoted as the "critical period" in the following, the measures are lifted on a gradually increasing rate, but initially (when the disease is not widespread in the population) only very slowly (R eff ≈ 1). It will be shown by analytical estimates, that the minimum duration for stable control during the critical period scales with C −1 0 , see Sec. 4.2. At the end of that period, the control function u opt (t) features a notable "wobble", see Fig. 4 (a, b) , where its curvature changes sign. Finally, the remaining measures are lifted on a gradually decreasing rate until herd immunity S (T ) /N < R −1 0 is established. Figure 4 (c) shows the trajectories of the controlled and the uncontrolled epidemic in different projections of the state space. By controlling the transmission of infection, the enormous excursion of the trajectory is prevented and the optimal path to a stable disease-free stationary state is taken. In particular, it should be noted that the stationary point where the uncontrolled epidemic stops is at S(T )/N R −1 0 (far in the stable regime), in contrast to the optimally controlled case S(T )/N R −1 0 , which is just sufficient to be stable. 

The results described above indicate that the state-dependent mortality rate f (C/C 0 ) effectively imposes a state-constraint that strictly enforces C < C 0 for P → ∞, i.e., a maximum number of simultaneously infected in a critical condition. In principle, this allows to investigate the optimal control for other (less extreme) scenarios, where e.g. the maximum number of simultaneously infected should be held far below the number of available ICUs (i.e., the meaning of C 0 will be reinterpreted). In this case, the increased mortality rate f 1 is an artificial parameter, that penalizes the excess of patients in a critical state above a freely chosen value of C 0 . By ramping up the control parameter P, an optimal solution with C (t) < C 0 for all t ∈ [0, T ] is found. The optimal control function computed in this way is completely independent of f 1 . Figure 5 shows the optimal control for different values of C 0 . The time course of the optimal effective reproduction number is qualitatively the same for all considered values of C 0 , see Fig. 5 (a, b) . Most notably, the time scale of the entire intervention scenario is governed by the duration of the critical period, during which the number of critical patients is held at C C 0 , see Fig. 5 (c). The full width half maximum (FWHM) time of that period T FWHM ∼ C −1 0 scales inversely with the peak number of simultaneously infected in a critical state, see Fig. 5 (d). The minimization of the disease-related deaths is controlled by the parameter P in the terminal cost function (5b). Figure 5 (e) displays the progression of the optimization routine into the targeted optimal state (i.e., without excess of C 0 ) while P is ramped up. At a certain value of P, which depends on C 0 , the routine reaches a plateau where both the number of disease-related deaths as well as the total costs of the intervention measures T 0 dt C (u (t)) become constant. The corresponding values of P, which correspond to the scenario that fully avoids excess of critically ill patients over C 0 , are located on that plateau and are marked by square symbols in Fig. 5 (e) . The optimized transmission function minimizes the number of disease-related deaths to a C 0 -independent value D min (T ) for P → ∞, but at total cost that scale with C −1 0 . The curve reflects the state-dependent fatality rate (2) . According to the present model, further reduction of disease-related deaths below D min (T ) can only be achieved by pharmaceutical interventions, in particular by vaccination. The result of the C 0 -independent number of deaths at the end of the epidemic is suspected to be an artifact of the simplified modeling framework, in which a macroscopic population with an averaged set of parameters is considered. Since mortality rates of individuals are typically strongly dependent on age and health condition, it might be advisable to extend the model and divide the compartments into several age or risk groups as in Refs. [1, 18] . The so-extended model features a matrix-valued transmission rate, which describes the infections caused by contacts within and between different groups, that could be further optimized by intra-and intergroup-specific selective measures. This is, however, beyond the scope of the current paper.

The numerically computed result illustrated in Fig. 4 (a, b) shows that during the critical period the populations S, R, and D change approximately linear, while all the others (active cases) are practically constant. To gain further insights, we consider the ansatz . The inset shows that the optimal control respects the stability requirement (12) during the critical period. The plot is for C 0 = 10 000.

(for t > t * )

where t * is a reference time that depends on the initial conditions, γ S , γ R , γ D are initially unknown rates and the infected sub-populations (E, I, H, C) ≈ (E * , I * , H * , C 0 ) are constant. From substituting the ansatz into the model equations (1), one obtains by a straightforward calculation analytical expressions for the rates

and the constants

Note that it holds γ S = γ R + γ D , i.e., the number of active cases remains constant since susceptible persons become infected at the same rate on which already infected either recover or die. The number of active cases in that dynamical equilibrium is found to be a multiple of C 0 :

Let us now come to the most important result of this section. The approach described above yields an instantaneous relationship between the current number of susceptible persons S(t) and the current optimal effective reproduction number, which is

for a certain range of t in t * < t < T crit with T crit defined below. This approximate relation is an interesting result, as it hints that the obtained optimal control steers the system's trajectory close to the stability boundary. Comparison with the stability criterion for the disease-free stationary state R 0 < N/S, see Eq. (16) , suggests that during the critical period one must make sure that

This allows to have a stable control of the number of actively infected cases, while the intervention measures can be gradually relaxed. Stable means that sufficiently small fluctuations of the number of infected are damped and do not lead to a new exponential outbreak of the epidemic. Indeed, substituting u (t) = (1 + ε) N/(R 0 S (t)) into the model equations (1) yields a linear, autonomous dynamical system (up to the state-dependent mortality rate (2)), which is easily seen to have a stable dynamical equilibrium for ε < 0. The understanding of the stability of the dynamical equilibrium during the critical period is the key insight to the qualitative understanding of the optimal control. The mechanism is described in more detail in Appendix B. The optimal effective reproduction number is plotted in Fig. 6 along with the analytical approximation (11) and the stability criterion (12) for the critical period. Because of the potentially enormous significance of this observation for the stable control of epidemics, this criterion shall be formulated once again in a different way: It is of utmost urgency to have the most accurate estimate of the cumulative number of cases at any point in time during the critical period. Since it holds S (t) = N − Ncases (t), where N is the total population and Ncases (t) is the cumulative number of cases that includes next to the active cases also the recovered and deceased population (12) can be reformulated as

Moreover, one can derive an upper limit for the rate of change of R eff (t) by differentiating Eq. (12) and using Eq. (1a) and the approximation I(t) ≈ I * (see above) aṡ

Hence, using the approximation (11), the rate on which the measures can be relaxed is limited by the square of the current value of the effective reproduction number R eff (t) times the rate γ S /N that is proportional to C 0 /N (see above), i.e., the share of the maximum allowed simultaneously critically ill patients on the total population. It holdṡ

The numerically calculated optimal control obeys the criteria (13)- (14) , see Fig. 6 , and is therefore (weakly) stable against smaller suddenly occurring new infections. The merely weak stability reflects the demand for minimal socio-economic costs, see Sec. 3. The two rules (13)- (14) for the optimal and stable steering of the effective reproduction number are expected to be generally valid, as they are widely independent of the details of the current model system. Equivalent results for a stable dynamical equilibrium with a constant number of infected cases are easily obtained for the much simpler SIR model, see Appendix B.

Finally, the characteristic duration T crit of the critical period can be estimated from Eq. (11) and the condition R eff (t * + T crit ) ≈ R 0 . One obtains

which is in excellent agreement with the numerically obtained values for the full width half maximum (FWHM) time plotted in Fig. 5 (d) .

Non-pharmaceutical measures to control the spread of infectious diseases and to prevent a potential collapse of the health care system must be precisely coordinated in terms of timing and intensity. Based on well-calibrated mathematical models, the optimal intervention strategy for specific scenarios and objectives can be computed using continuous-time optimal control theory. In this paper, an extended SEIR model was calibrated to reproduce the data of the initial exponential growth phase of the COVID-19 epidemic in Germany. Optimal control theory has been applied for the scenario in which an effective vaccine is impossible or will never be found and the epidemic must be controlled with purely non-pharmaceutical measures. The computed optimal transmission rate describes a single-intervention scenario that satisfies competing constraints: First, the minimization of the diseaserelated deaths by strictly avoiding an overflow of intensive care resources and, second, the suppression of a second epidemic wave by establishing sufficient natural immunity at the end of the measures. Moreover, the total costs of the intervention were kept at a necessary minimum for socio-economic reasons.

The optimal control computed in this paper is based on a relatively simple epidemiological model and the quantitative results for the COVID-19 epidemic depend to some extent on the chosen parameters. Therefore, the quantitative results must be treated with caution and should be compared with the results from more realistic models. Nevertheless, the described intervention scenario was shown to be qualitatively robust under parameter variation and has well-understood stability properties. The comparison of the optimal solution with the stability criteria (13)-(14) reveals, however, that the obtained solution is in fact a tightrope walk close to the stability boundary of the system, where socio-economic costs and the risk of a new outbreak must be constantly balanced against one another. Furthermore, our analysis clearly shows that the goal of achieving herd immunity via natural infections is either extremely expensive (in terms of socio-economic costs due to measures maintained over a long period of time) or extremely dangerous (due to the constantly high load on intensive care resources just below the stability limit). Note that the values of C 0 considered in the computations are very high throughout. In any case, in view of the long duration and the enormous number of infections that this route entails, as well as the uncertain role of sequelae and the uncertain prospects for appropriate vaccines, it is strongly advisable to consider other strategies, in particular the attempt to eradicate the epidemic completely [33] . Fig. 7 The stability of the disease-free stationary state depends on the size of the susceptible populationS and the basic reproduction number R 0 . ForS/N < R −1 0 , the outbreak of an epidemic is suppressed by a sufficiently high degree of herd immunity.

Without intervention, i.e. u (t) = 1, the system (1) has a family of disease-free stationary states x = S , 0, 0, 0, 0,R,D . The stability of a stationary state with respect to small perturbations x →x + δx (t) is determined by the sign of the real parts of the eigenvalues η of the linearized system's coefficient matrix

From the characteristic polynomial

one obtains the eigenvalues

and the threefold degenerate eigenvalue η (0) = 0. Since c (1 − f 0 ) < 1, it holds η

± < 0. The leading eigenvalue is η (1) + , which is negative for

see Fig. 7 . Hence, the disease-free stationary state is unstable if the susceptible population size exceeds a critical threshold value that is given by the inverse basic reproduction number (3) . ForS/N < R −1 0 an epidemic outbreak is suppressed by a sufficiently high degree of herd immunity.

This section sheds more light on the key mechanism behind the control during the critical period, see Sec. 4.2, where the number of active cases is kept constant while the intervention measures are gradually relaxed. In order to highlight the supposed general validity of the rules (13)- (14) , here the much simpler SIR (susceptible-infectious-recovered) model is considered. The model equations read [3] 

For now, this is a nonlinear and non-autonomous dynamical system. Our goal is to keep the number of simultaneously infected active cases constant, i.e., we seek for a solution withİ = 0. Then, from Eq. (17b), we deduce

where R 0 = β/γ i and R eff (t) = R 0 u (t). Substituting this into the SIR system's model equations (17) (18), which describes the strength of the non-pharmaceutical measures, we obtain

Initially, the effective reproduction number is R eff (0) = 1 and for t > 0 the value is monotonously growing. The intervention stops at the time T , where all measures are lifted an the effective reproduction number equals the basic reproduction number R eff (T ) = R 0 . The rate of change of the effective reproduction number is obtained aṡ

i.e. it grows with increasing speed, depending quadratically on the current value of R eff (t) and the rate γ S /N ∝ I * /N , which is proportional to the (constant) share of the infected cases on the total population. Finally, we study the stability by investigating the impact of a small shift of the control function

This yields a linear system in the vicinity of the dynamical equilibrium:

The solution I (t) = I (0) exp (εγ i t) shows that for ε > 0 the number of infected grows exponentially fast, while for ε < 0 it is exponentially damped. Thus, for ε < 0, i.e., for

which is the criterion (12), the control scenario is robust with respect to small perturbations of the number of cumulative cases, e.g., small spontaneously occurring infection foci or small numbers of imported cases. The second criterion (14) for stable control follows immediately from Eq. (12) and Eq. (17a) or (1a), respectively, as given in the main text. 10 -1 number per disease-related deaths Fig. 8 Ratio of the state variables in the initial exponential growth phase and the number of disease related deaths. The numerically exact solution (solid lines) is plotted along with an analytical approximation (solid lines) that holds in the early stage of the epidemic. The corresponding algebraic relations are used to describe ratios between different sub-populations to facilitate the parameter adjustment. Symbols indicate the reported number of diseaserelated deaths (black), estimated number of cases (grey) and estimated ICU load (purple) of the COVID-19 epidemic in Germany. See the text for details. During mid-March, strict social distancing measures were implemented, that flattened the initial exponential growth.

The parameters are adjusted such that the model reproduces the data of the early exponential growth phase of the COVID-19 epidemic in Germany. It is of course questionable to calibrate an epidemic model to a single country, but in a scenario with extensive border closures this seems to be justified. In the exponential growth phase of the epidemic, all sub-populations grow exponentially with the same rate, see Fig. 2 (b). This observation can be exploited to derive a series of algebraic equations (which hold approximately in the initial phase of the epidemic) that relate all state variables to each other. On the basis of empirical data (reported number of cases and deaths etc.), several missing model parameters can be directly determined from the algebraic relations. The number of reported cases and deaths used in this study is based on the figures provided by the Robert Koch-Institute [20, 24] .

One starts with the ansatz

where Γ is the initial exponential growth rate that is estimated from reported data (see Fig. 2 

and the relation between the growth rate and R 0 :

Note that Eq. (20) is equivalent to the equation for the leading eigenvalue η (1) + if the whole population is susceptible, i.e. Γ = η (1) + |S =N (see Appendix A). Hence, Eq. (20) implies that the exponential growth rate Γ changes sign at R 0 = 1, i.e., the epidemic recedes for R 0 < 1. The mean incubation period was reported to be 5.1 d, but there are indications that the latency time may be shorter [16] . Assuming the onset of infectiousness 2.5 d before the onset of symptoms, this implies an average latency period of γ −1 l = 2.6 d, i.e., the latency period is assumed to equal roughly half of the incubation period. The reported values of the basic reproduction number R 0 are heavily scattered. According to the Robert Koch Institute, serious estimates range between 2.4 and 3.3 [25] . In the following R 0 = 2.7 shall be used, which is situated approximately in the middle of the interval in question. From Eq. (20), the corresponding average infectious period is obtained as γ −1 i ≈ 2.35 d. The overall infection fatality rate of COVID-19 was estimated as 0.66 % [32] , such that (1 − m) cf 0 = 0.0066. On April 8, the Robert Koch Institute reported that a fraction of f 0 = 0.31 patients in a critical state died (without ICU overflow) [26] . Finally, the fraction of infected with at most mild symptoms is estimated as m = 0.92, such that c = 0.0213/ (1 − m) ≈ 0.266.

Substituting the exponential ansatz (19) in Eqs. (1d)-(1g), yields

The analytically obtained ratio between all subpopulation and deaths (which are believed to be the most reliably reported data) are plotted along with the corresponding numerically exact result for the initial uncontrolled epidemic in Fig. 8 (b) . The analytical results imply the relation

Unfortunately, there is only little data available on the demand for ICUs in the early phase of the epidemic. In mid-March 2020, i.e. near the end of the initial exponential growth phase, the German Interdisciplinary Association for Intensive Care and Emergency Medicine (DIVI) initiated a register that reports on the availability of ICUs in Germany [10] . On March 27, 687 out of 1 160 hospitals with ICUs contributed to the register and reported a total number of 939 COVID-19 patients in a critical state receiving intensive care [27] . At the same day, 253 disease-related deaths were reported. From the estimated ratio C/D ≈ 6.3 (the actual number of critical patients was estimated based on the ratio of contributing and non-contributing hospitals as C ≈ 1 586), the average period after which patients in a critical state either recover or die, is estimated from Eq. (21) as γ −1 c ≈ 7.5 d. Finally, assuming that only r = 2/3 of all cases have been discovered initially and an assumed average time delay between infection and report of cases of ∆tr = 5 d, the number of actual cases is estimated from the number of reported cases as N est cases (t) = r −1 N rep cases (t + ∆tr) = r −1 e Γ ∆tr N rep cases (t) ≈ 5.5 N rep cases (t). This yields a good agreement between the simulated number of cases (Ncases = E + I + H + C + R + D) and N est cases before measures came into force, see Fig. 2 (a) and Fig. 8 . The average time between infection and death ∆t d can be estimated from the ratio N est cases (t) /D (t) ≈ 2370 (see Fig. 7 (b)) and N est cases (t − ∆t d ) = N est cases (t) e −Γ ∆t d = D (t) as ∆t d ≈ 29.9 d. The equations for R and D are decoupled and can be solved via direct integration once the solution of the remaining system is known. The corresponding Lagrange multipliers are λ R (t) = 0 and λ D (t) = P. The initial time conditions guarantee u (0) = 1 (no intervention, see Eq. (8)) at the beginning of the scenario.

A first study on the impact of current and future control measures on the spread of COVID-19 in Germany

A hybrid multi-scale model of COVID-19 transmission dynamics to assess the potential of non-pharmaceutical interventions

Compartmental models in epidemiology

Modelling transmission and control of the COVID-19 pandemic in Australia

Optimal COVID-19 epidemic control until vaccine deployment. medRxiv

Projecting the spread of COVID19 for Germany

Sequential data assimilation of the stochastic SEIR epidemic model for regional COVID-19 dynamics

Impact of non-pharmaceutical interventions (NPIs) to reduce COVID-19 mortality and healthcare demand

Estimating the number of infections and the impact of non-pharmaceutical interventions on COVID-19 in 11 European countries. Imperial College London

German Interdisciplinary Society for Intensive Care Medicine (DIVI): DIVI Intensivregister

The mathematics of infectious diseases

Modeling the control of COVID-19: Impact of policy interventions and meteorological factors

A contribution to the mathematical theory of epidemics

Estimate of the development of the epidemic reproduction number Rt from coronavirus SARS-CoV-2 case data and implications for political measures based on prognostics

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

The Incubation Period of Coronavirus Disease 2019 (COVID-19) From Publicly Reported Confirmed Cases: Estimation and Application

Optimal control

COVID-19: Development of a robust mathematical model and simulation package with consideration for ageing population and time delay for control action and resusceptibility

Optimal timing for social distancing during an epidemic. medRxiv (2020)

Intensive care management of coronavirus disease 2019 (COVID-19): Challenges and recommendations

Coronavirus: The hammer and the dance

SARS-CoV-2 Steckbrief zur Coronavirus-Krankheit-2019 (COVID-19)

Täglicher Lagebericht des RKI zur Coronavirus-Krankheit-2019 (COVID-19)

Täglicher Lagebericht des RKI zur Coronavirus-Krankheit-2019 (COVID-19)

Adjoint-based data assimilation of an epidemiology model for the Covid-19 pandemic in 2020

Solving ODEs with Matlab

Flattening the curves: on-off lock-down strategies for COVID-19 with an application to Brazil

Modeling, state estimation, and optimal control for the US COVID-19 outbreak

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

Stellungnahme der Helmholtz-Initiative "Systemische Epidemiologische Analyse der COVID-19-Epidemie

Modelling the potential health impact of the COVID-19 pandemic on a hypothetical european country

A new coronavirus associated with human respiratory disease in China

A novel coronavirus from patients with pneumonia in China