Intensive care unit occupancy predictions in the COVID-19 pandemic based on age-structured modelling and differential flatness

The COVID-19 pandemic confronts governments and their health systems with great challenges for disease management. In many countries, hospitalization and in particular ICU occupancy is the primary measure for policy makers to decide on possible non-pharmaceutical interventions. In this paper a combined methodology for the prediction of COVID-19 case numbers, case-specific hospitalization and ICU admission rates as well as hospital and ICU occupancies is proposed. To this end, we employ differential flatness to provide estimates of the states of an epidemiological compartmental model and estimates of the unknown exogenous inputs driving its nonlinear dynamics. A main advantage of this method is that it requires the reported infection cases as the only data source. As vaccination rates and case-specific ICU rates are both strongly age-dependent, specifically an age-structured compartmental model is proposed to estimate and predict the spread of the epidemic across different age groups. By utilizing these predictions, case-specific hospitalization and case-specific ICU rates are subsequently estimated using deconvolution techniques. In an analysis of various countries we demonstrate how the methodology is able to produce real-time state estimates and hospital/ICU occupancy predictions for several weeks thus providing a sound basis for policy makers.

groups. By utilizing these predictions, case-specific hospitalization and case-specific ICU rates are subsequently estimated using deconvolution techniques. In an analysis of various countries we demonstrate how the methodology is able to produce real-time state estimates and hospital/ICU occupancy predictions for several weeks thus providing a sound basis for policy makers.

I
Reported active cases [reported number of currently infected people (positive tested)] S Susceptibles (people that are susceptible to the COVID-19 disease) u Aggregated exogenous drivers (leading indicator for the course of the disease obtained from nonlinear control theory) N Population size (number of people in the considered country or age group) β Transmission rate (rate (speed) describing the spread of COVID-19) γ Recovery rate (rate at which infected individuals recover) π Production rate (daily fow of individuals from compartment S to I, i.e. π = βSI/N for the standard SIR-model) r Case-specific admission rate (share of people tested positive with COVID-19 who are admitted to the hospital or ICU) χ Occupancy (number of people/occupied beds in hospital or ICU at a specific day) d max Maximal days in hospital or ICU (maximal number of days an individual of an age group occupies a hospital or ICU bed due to COVID-19) ϕ(τ ) Probability of care (probability for being in care τ days after infection, conditional to prior hospitalization or ICU assignment) Indices a Age group (placeholder for age groups) h Hospital (refers to hospital) icu Intensive care unit (refers to ICU) 1 Younger age group (refers to the younger age group of a country) 2 Older age group (refers to the older age group of a country)

Introduction
The ongoing fight against recurring epidemic waves of the SARS-CoV-2 pandemic is a complex undertaking that requires several key directions of action such as: (i) social distancing and personal protective equipment, (ii) testing for infection, (iii) quarantine of infected people and (iv) vaccination. Management of healthcare systems which are treating the actual COVID-19 patients remains a central issue during the ongoing pandemic. The necessity to avoid overloading of the healthcare system is imperative as the level and quality of medical care and particularly mortality are directly related to the available capacities. Within the healthcare system, the ICU capacity constitutes a major bottleneck for policymakers to take appropriate measures.
Therefore, predictions of admissions to hospitals and intensive care units (ICUs) are a critical factor in the decision-making process. They are equally important in the phases of upcoming surges and in phases of decline, where far-reaching measures such as lockdowns should be lifted as quickly as possible.

Occupancy models
The dynamics of the epidemic along with the time lag between the infection and the actual hospitalization as well as the duration of stay (ranging from a few days to several weeks) make it difficult to predict the occupancy of COVID-19 patients in hospitals and ICUs in particular.
To infer future hospital and ICU occupancies from reported case numbers, in this paper an approach that comprises the online estimation of case-specific hospitalization/ ICU admission rates is proposed using deconvolution techniques. It explicitly accounts for time lags between infection and hospitalization as well as for the distribution of the length of stay. To allow for this approach to make predictions depending on different future evolutions of infections, an age-structured compartmental model is proposed along with a methodology that provides estimates of the states of that model and the otherwise unknown exogenous inputs driving its dynamics requiring only the reported infection data.

Modelling approach and review
The compartmental model utilized in this work is based on the classical SIR-model [1] but has been extended recently by an additional compartment ('contact-less') in combination with an approach to determine unknown exogenous input u(t) driving its dynamics [2]. The latter is based on the concept of differential flatness [3,4]. The exogenous input u(t) aggregates all unknown drivers of the epidemic and thus enables the realistic analysis and prediction of epidemiological dynamics, in particular recurring waves.
The methodology in [2] is now further augmented by segregation into two coupled age groups in order to account for the fact that both hospitalization and ICU occupancy are strongly age-dependent. Age-structured compartmental models are commonly used in epidemiology. For example, the recent publication [5] proposes two age groups to model the spread of tuberculosis. A more complex stochastic compartmental model of the COVID-19 pandemic is utilized in [6] where five age groups are employed. Utilizing data for Belgium they also predict hospitalization and ICU occupancy. No estimation of external drivers of the epidemic is given, though, and the prediction accuracy for new hospitalizations can hardly be determined based on the pre-sented validation. Data from the COVID-19 pandemic in Switzerland is treated in [7]. A discrete model based on compartments with explicit duration of the individual phases is presented. Seventeen age groups are distinguished, and both hospitalization and ICU occupancy for these groups are computed. The model is parametrized with the measured age distribution in the hospitals. Although several simulations are shown to analyze the effects of different measures, no methodology for real-time predictions is presented. A compartmental model with four age groups is presented in [8] which focuses on an optimal control approach and presents simulations of closed-loop scenarios.
The topology of the population is also an important factor to impact the epidemics. In [9] an extended SIR model is presented which describes the investigation of the epidemic spread in complex networks by considering the propagation vector and infection delays. In this context, [10] studies the co-evolution of multiple information dissemination and epidemic under the influence of mass media.
A Monte-Carlo SEIR model is employed in [11] to model and simulate COVID-19 outbreaks in the Netherlands, South Korea, China and Italy. Although the model is used for forecasts of case numbers and hospitalizations, ICU occupancy is not treated but the forecasts could not be validated. In [12] a SEIAR model is presented with four distinct age groups for the COVID-19 pandemic in Hunan and Jilin provinces. The model was parameterized based on literature and own data fitting, respectively. The main result is the estimation of SAR-values (secondary attack-rate values) between different age groups.
In [13] a methodology to forecast the patient numbers admitted to hospital and ICU is presented. The model is restricted to data from the hospital and no generic pandemic model is included, but time-varying admission rates are estimated. Forecasts up to five days as well as the respective validation data are presented. In [14] the hospitalization rate is estimated for New Mexico by pre-filtering and dividing reported hospitalization cases by reported COVID-19 case numbers; no pandemic model is employed. Several different models are compared and analyzed in [15] for forecasting ICU occupancy. The model outputs are merged by a trimmed-mean approach to predict the occupancy directly. Predictions and validation data is presented for two-week forecasts, however, age segregation is not considered. Another model without age segrega-tion is presented in [16]. ICU admission rates, time lag until ICU admission and duration of intensive care are directly estimated. Predictions show squared correlation values between 0.00 and 0.99 depending on the data set, the better values being achieved during exponential growth phases. Data from a publicly available forecasting tool [17] are processed to obtain statistical estimates of occupancy. This statistical model is then trained to minimize the prediction error. No age segregation is considered and especially during peak values of the pandemic model deviations are considerable. The trimmed mean of autoregressive, machine learning, and epidemiological models is utilized in [15] for estimating the ICU occupancy. One-and two-week predictions are shown together with reported cases, but no distinction is made between age groups.
A different approach is shown in [18], using static and dynamic incidence models as well as a care pathway model to generate predictions for ICU and hospital beds. Using local data from health districts in Germany, forecasts were made several weeks into the future. However, these were compared to the actual cases only for the first 23 days, mostly falling outside the 95% confidence interval. A differentiation of the hospital data into age groups was also not considered here. A prediction of hospital and ICU beds was also performed by [19], using an adjusted SEIRD model for its determination. A two-week forecast was shown for the Mexico City metropolitan area, which was also compared to actual cases. However, no predictions were shown at other time periods, making it difficult to draw conclusions about the performance of the model. In addition, no age classification was considered.
Besides the hospitalization numbers, the casespecific hospitalization and case-specific ICU admission rates r H and r ICU , respectively, constitute decisive latent variables in the analysis of the epidemic. In the method presented in this paper both case-specific admission rates are treated as variable over time and a methodology is proposed to estimate them in real-time using deconvolution techniques. To this end, only the reported active cases along with hospital and ICU occupancy data are required. Another approach to model ICU occupancy is presented in [16]. However, the casespecific admission rates are not estimated by deconvolution but are directly taken from literature or estimated externally.  Fig. 1 Schematics of the proposed methodology to estimate and predict the occupanciesχ for hospital and ICU

Methodology
In the methodology proposed in this work, the predictions of hospital and ICU occupancy, respectively, are based on the predictions of active COVID-19 case numbers using the flatness approach [2] along with projections of the case-specific admission rates, cf. Fig. 1. Predictions of hospital and ICU occupanciesχ h andχ icu are obtained from a forward dynamic occupancy model, which by itself uses predictions of the states of an forward epidemic model (i.e. the number of infectedÎ , susceptiblesŜ) and case-specific admission ratesr h,icu .
To this end, in Sect. 2, an age-structured compartmental epidemic model is introduced and it is shown how the unknown exogenous inputs u that drive its dynamics can be estimated. Subsequently, it is demonstrated how these inputs are utilised to obtain predictions of the future epidemiological statesÎ (t + τ ) and S(t + τ ) using only information available up until time t. These predictions are then used for the occupancy predictionsχ(t + τ ).
In the proposed approach the case-specific hospitalization and ICU admission ratesr h andr icu are treated as unknown latent variables. Besides the infection dynamics itself, they are the second major factor that seriously impacts the burden on ICUs and hospitals. In Sect. 3 it is described how these admission rates can be estimated in real-time based on deconvolution techniques.
Finally, dynamic occupancy forecastsχ(t + τ ) for different countries are presented and discussed in Sect. 4.

Age-structured compartmental model with exogenous drivers
The description of epidemics is often based on compartmental models [20] which are versatile approaches and can easily be extended or adapted, as shown in [21][22][23]. In this section, the compartmental model with exogenous drivers as presented in [2] is recaptured and then augmented to a discrete, age-structured model. The estimates of the epidemiological states obtained from the proposed segregated modelling approach provide the basis for estimating the hospital and ICU occupancy.

Compartmental models with exogenous drivers
Compartmental epidemiological models constitute a set of coupled autonomous ordinary differential equations. The dynamics are exclusively determined by the initial states and the parameters. However, the COVID-19 pandemic is significantly driven by exogenous inputs. There are different approaches to extend compartmental models to non-autonomous systems with exogenous inputs, e.g. [24][25][26], so that they can accurately describe multiple epidemic waves. In this work, we focus on the methodology introduced in [2] which is based on the mathematical property of differential flatness. Various compartmental models have the property of differential flatness [27,28]. Extending them with exogenous drivers u(t) [2] is briefly discussed for the SIR model: with I as the number of infected or active cases, S as the susceptibles, R as the recovered individuals and N as the population size. The transmission rate β and the recovery rate γ are parameters of the system. The exogenous driver u(t) can be obtained by differentiating system (1) with respect to time. Due to the extension of system (1) with exogenous drivers, the assumption that N is the sum of compartments S, I, R does not hold. Therefore, the additional compartment of contact-less C is introduced. The population size is thus where C has no additional state equation. Differential flatness entails that u(t) and S(t) can be expressed solely by I (t) and a finite number of its derivativesİ andÏ , as will be shown for an age segregated CSIR model in the sequel.

Age-structured CSIR compartmental model
In order to account for how much more severely COVID-19 affects the elderly, it is appropriate to segregate the CSIR model into two discrete but coupled age groups, indicated by "1" and "2": In Equations (2-7), the compartments S, I and R are split into the two age groups. The age limit for the segregation into these groups is 65 years for Austria but can vary depending on the investigated country and the available data. Equations (2-7) are visualised in Fig. 2, where the flows of individuals between the compartments and their couplings, including compartments C 1 and C 2 , are shown. Notably, each age group is driven by its own exogenous input u 1,2 making it a multi-input system.
The system represented by (2-7) utilises the three transmission rates β, β 12 and β 21 which describe transmissions within and between the age groups. The transmission rates and the recovery rate γ constitute the parameters of the compartmental model.
For the estimation of the ICU and hospital occupancy presented in the following sections, age group specific production rates are introduced: Here, λ 1,2 represent the forces of infection for the specific age groups and S 1,2 the susceptibles, respectively. Based on equations (2-5), they are given by which establishes a cross-coupling of the dynamics of the two age groups.

Estimation of unknown exogenous drivers
Much like the basic SIR model (1) the age segregated model has the property of differential flatness, differing only in the aspect that now a multi-input multi-output system is considered [29,30]. Specifically, using the reported case numbers I 1 (t) and I 2 (t) and their first and second derivatives, respectively, one can derive algebraic equations for the number of susceptibles S 1,2 : The exogenous drivers u 1 and u 2 are obtained by deriving equations (4-5) with respect to time. After replac-ingṠ 1,2 with equations (2-3) and S 1,2 with (10), the exogenous drivers u 1,2 for the discrete age groups are obtained as: inverse epidemic model forward epidemic model . 3 Top: Using the differential flatness property, an inverse use of the epidemic model provides real-time estimates of its states as well as its exogenous inputs that drive the dynamics. Bottom: Using projections of the exogenous inputs, the compartmental epidemic model provides predictions of its future states Using the time course of u 1 (t) and u 2 (t) determined from (11), one could reproduce the observed dynamic behaviour of the epidemic for each age group exactly. In this sense, (11) can be interpreted as the inverse of the epidemic model, Fig. 3 (top). How u 1 (t) and u 2 (t) can be used to describe and forecast the pandemic is addressed in the following section.

Analysis and forecasts based on exogenous drivers
A real-time analysis of the COVID-19 epidemic of a country involves the determination of its missing states S 1,2 (t) as well as its exogenous drivers u 1,2 (t) using the equations outlined in the previous section. Conceptually this can be interpreted as an inversion of the age segregated inverse epidemic model, Fig. 3 (top).
Once the current states of the model are available, its future state trajectory is solely defined by the future course of its exogenous drivers (e.g.û(t + τ )). Hence, forecasts of the epidemic can be obtained by suitable projections ofû 1,2 (t + τ )) which are then fed into the epidemic model. This amounts to a forward use of the model, Fig. 3 (bottom). The projectionsû(t + τ ) therefore constitute the inputs to predictÎ (t +τ ) andŜ(t +τ ) based on Equations (2)(3)(4)(5).
Obviously, many different methods could be considered how to obtain the necessary projections ofû(t +τ ), e.g. regression models or expert knowledge. One suitable approach that resulted in sufficiently accurate projections throughout the pandemic is to apply gradient based linear extrapolations forû. The strength of this relatively simple approach is its straightforward adaptability to already observable trends in u(t) and they proved to be a robust and especially transparent method. This is illustrated in Fig. 4 by means of rolling three-week forecasts (i.e. 0 < τ ≤ 21) for the active case numbers for Austria. The first age group corresponds to the younger population, aged below 65 years, whereas the second group is related to the elderly at or above 65.
The first two subplots of Fig. 4 show the reported numbers of infected for the two age groups I 1,2 (solid lines). These are to be compared against the rolling forecastsÎ 1,2 (t + τ ) (dashed lines) for different points in time. The beginning of each forecast is marked by a circle. The estimated aggregated exogenous drivers u 1,2 (t) obtained from the inverse epidemic model are displayed as solid lines at the bottom of Fig. 4, together with the linear projections ofû 1,2 (t +τ ) (dashed) which are required to obtain the case number forecasts. It can be observed, that the actual courses of the infected are predicted accurately throughout each three-week prediction phase for both age groups in most cases. Notably, a correct prognosis in the vicinity of the peak in fall 2020 in Austria proves to be difficult as the prediction range contained an unpredicted lockdown. The corresponding results for France and the Netherlands are given in the Appendix in Figs. 11 and 15, respectively.

Estimation of case-specific hospitalization and case-specific ICU admission rates for age groups
In addition to the number of infected people, a key indicator for the severity of the epidemic situation of a country are the remaining capacities of hospital and ICU beds. The methodology to predict the occupancies is based on the case-specific hospitalization/ ICU admission rates. They are defined as the share  [31] of patients who are e.g. admitted to ICU following an infection. While in many related works case-specific admission rates are treated as given and constant over time [16,32,33], we propose an approach where these rates are time variant and hence have to be determined from the available data. Especially in view of the upcoming phase of COVID-19 where vaccination rates in many countries are on a steady increase, case-specific admission rates serve as important latent variables that can give valuable hints about the success of vaccination campaigns.

Hospitalization and ICU occupancy model
An important aspect for modelling and estimating the occupancies is the available data. Particularly important is data that links the length of stay durations in hospital and ICU care to the day of admission. For agestructured models, this information needs to be reported in an even higher resolution. If these data are available, convolution based approaches to model the occupancies can be considered [16].
To this end we define ϕ(t − θ) as the conditional probability of being in care on a certain day t, conditional to prior admission a number of days (t−θ ) earlier.
Using the available hospitalization statistics, ϕ(t − θ) can be determined in a straightforward way. Similar to [16], convolution techniques are used to link the current production rate of infections π(t) in (8), the conditional probabilities ϕ(t − θ) as well as the associated casespecific admission rates r (t) to determine occupancies.
The approach to determine the occupancies χ a for the two age groups a=1,2 resulting from the above considerations is expressed by the discrete convolution for hospital beds ("h") and intensive care units ("icu"), respectively. In (12), π a is the number of new infection cases of the respective age group according to (8), r h,a and r icu,a are the time varying case-specific hospi-inverse dynamic occupancy model (deconvolution) . 5 The case and age specific admission ratesr a are estimated in real-time using a dynamic occupancy model in combination with deconvolution techniques. Therein, the states of the epidemic model (i.e. I a , S a ) and the reported occupancies χ a are utilized talization and case-specific ICU rates. The maximum number of days individuals are considered to be in care is denoted by d max .

Estimation of case-specific rates based on deconvolution
In (12), r a can be seen as a latent variable which has an important meaning for the pandemic: As will be shown later (cf. Fig. 6) the case-specific admission rates show significant fluctuations over time. This can be due to e.g. vaccination effects or new trends in medical treatment policies. In this work, the case-specific admission rates are therefore treated as time-varying and are estimated based on available measurements as presented schematically in Fig. 5. Much like the exogenous inputs that drive the epidemic, estimation of the case-specific admission rates is based on the principle of dynamic inversion, whereas this time specifically the occupancy model is considered. Based on the epidemiological states I a (t), S a (t) from (2-5) and the reported occupancies (e.g. χ h,a (t)) estimates of the case-specific admission ratesr a (t) are obtained using deconvolution techniques. These are also known under the names polynomial division or backsolving [34][35][36]. For this matter, a Toeplitz matrix is introduced which contains the conditional probabilities ϕ(t) as for each respective age group 1 . Second, an n × 1 vector p containing the production rates π(t j ) at each observation time t j ( j = 1, . . . , n) is defined as The last observation at t n might be the day the estimation of r itself is conducted, if data are available up to this day. In the same fashion the n × 1 vector of observed occupancies and the yet unknown vector of case-specific admission rates r Using the above definitions, under ideal conditions the observed occupancies are formally related to the respective case-specific admission rates through which is the equivalent of (12), evaluated for all observed time instants and expressed in vector-matrix notation. There are, of course, several ways to determine an estimater of the vector of case-specific admission rates, whereby each of them is inherently linked to the minimization of some kind of error criterion, like, e.g.
When practically computingr from reported data it turns out that for several countries the result can be sensitive to fluctuations which stem, e.g. from irregularities in the data reporting process. To overcome this problem, regularisation of the second time derivative of r (t) proved to be a suitable remedy. The motive for this lies in the nature of the time evolution of admission Fig. 6 Case-specific hospitalization and ICU admission rates in Austria for both age groups (solid lines), including rolling forecasts (dashed lines). A case-specific admission rate of 100% means, that every infected person has been hospitalized as well. The forecastsr start on the same days as the projections ofû in Fig. 4 rates: Even though they are time variant their fluctuation over time and in particular their second derivatives can be associated with, e.g. large scale changes in treatment policies applied in hospitals or with the severity of the cases. As each of these underlying trends is subject to certain smoothness, it is reasonable to introduce regularisation. Then,r is obtained from where the matrix H is chosen such that the regularisation term penalises a measure corresponding to the second time derivative ofr. A proper choice is Finally, Q is a diagonal weighting matrix that controls the amount of regularisation.
The solution of (16) can be obtained analytically, yielding Equation (17) is evaluated separately for each age group as well as for the hospital and ICU occupancy, respectively. This results in four estimatesr * h,1 ,r * h,2 , r * icu,1 andr * icu,2 for each country. The estimated case-specific admission rates obtained from (17) are presented for Austria in Fig. 6.
They are shown as solid lines in all plots. Black lines indicate the younger age group (0-64 yrs.) whereas red lines correspond to the older group (65+ yrs.). It can be seen, that the case-specific hospitalization and ICU admission rates for all age groups vary significantly over time. It can be also clearly seen that the risk of being assigned to hospital or ICU care after an infection differs by roughly one order of magnitude between the age groups. The corresponding results for France and the Netherlands are again given in the Appendix.

Analysis and forecasts of age-structured hospital and ICU occupancies
The following shows how the occupancies of hospital and ICU beds for different age groups based on the estimated exogenous drivers and the estimated casespecific admission rates can be predicted and how rolling forecasts are obtained. The graphs in Fig. 6 show the projected case-specific admission ratesr (t + τ ) in addition to the estimated admission ratesr * (t). The beginnings of these projections are marked with an orange circle and all projections last three three weeks into the future (i.e. 0 < τ ≤ 21) making use of extrapolation by low-order polynomials. Since the projections ofr (t + τ ) and the epidemic states are required for the occupancy forecasts, the same points in time are chosen as for the projections ofû in Fig. 4. Figure 7 schematically shows how the predictions of hospitalχ h,a (t + τ ) or ICU occupanciesχ icu,a (t + τ ) are determined.

Hospitalization and ICU occupancy predictions
Predictions about future hospitalizations and ICU occupancies are obtained from the convolution sums (12) where the upper limit is now set to the respective future day t + τ :

Fig. 7
Using the convolution sums in (12), predictions of occupancies are obtained from predicted states of the epidemic (i.e. I a ,Ŝ a ) and projections of case specific admission ratesr a Note that the convolution sums in (18) have to be segregated into two parts each: The first part, which contains only those time arguments θ which refer to past time instants up to the present (i.e. t) and the second part, where only future time instants are considered. Consequently, the first sum, e.g. for hospitalizations, uses r h,a and π a , whereas in the second sum the respective predictionsr h,a ,π a are used. The latter predictions of the production rates (8) are given aŝ which in turn make use of the predictions of the epidemic statesÎ a (θ ),Ŝ a (θ ), see Fig. 3.

Analysis of hospitalization and ICU occupancy forecasts
In the case of Austria, the predicted hospital and ICU occupancies for both age groups obtained from (18) are shown in Fig. 8, for France and the Netherlands the reader is referred to the Appendix. In the corresponding Figures, the three week predictions for τ = 1, 2, . . . , 21 days ahead are represented by dashed lines in the colour of the respective age group (black, red). For comparison, the actually reported numbers are shown as solid curves. Note that the starting points of the predictions, displayed as orange circles, were again chosen equal to those in Fig. 4 in order to keep context. Naturally, prediction quality is of especial interest during strong epidemic surges as the one that happened in Austria during fall 2020. To shed more light on this particular phase, a more detailed analysis is shown in Fig. 9, where predictions were updated every four days. Following the assumption that the epidemiological course will significantly change with altering governmental interventions, predictions extending into such a turn of events (i.e. lockdowns) are cut off.

Statistics of occupancy predictions
This subsection presents an analysis and statistics for multi-week forecasts of the age-structured hospital and ICU occupancies for Austria. Corresponding results for France and the Netherlands can be found in the appendix. Figure 10 shows the relative errors (with respect to the corresponding highest reported historical occu-   To allow for a statistically meaningful evaluation, predictions were generated/updated on a daily basis from August 1, 2020 to May 11, 2021. In this analysis, two different variants are considered: First, the predicted occupanciesχ h,a (t + τ ) and χ icu,a (t + τ ) defined by (18) are analysed. In a second variant, in comparison,π a (θ ) in the second sum of (18) was replaced by π a (θ ) to examine the impact of epidemic forecasts: Obviously, the sums in (20) can only be evaluated expost, once the epidemic data up until t +τ are available. Although this variant cannot be used for forecasts in real-time, the retrospective evaluation ofχ h,a (t + τ ) andχ icu,a (t + τ ) based on the observed production rate enables the accuracy of the dynamic occupancy model and the admission rate prediction to be quantified. In general, the comparison of the two variants allows to investigate the sensitivity of the predictions to changing case-specific admission rates over the predicted time interval as well as the influence of the projections of future epidemiological development based on the exogenous drivers.
The representation of the prediction errors in normalised histograms (orange bins forχ h,a andχ icu,a , blue bins forχ h,a andχ icu,a ) in Fig. 10 allows for an analysis of their distribution and provides a means to determine 95% confidence intervals that can be used to assess the accuracy and reliability of future rolling forecasts. While the relative prediction errors after 7 days are small for both, more significant deviations with larger errors for the (true) predictions based on (18) appear after 14 and 21 days, as expected. The higher  (20), i.e. the prediction under the assumption that the future course of the infection numbers is already known, shows that the dynamic occupancy model in combination with the admission rate prediction provides very accurate forecasts over several weeks.
The results of the statistical analysis are also presented in Table 1 in terms of the RMSE and in Table 2 using the normalised RMSE (NRMSE), respectively. In particular, the NRMSE values (again, related to the highest observed historical occupancies) of the prediction errors show that the proposed methodology can produce very accurate predictions of future hospital and ICU occupancies.

Conclusion
We presented a methodology to quantitatively predict age segregated COVID-19 case numbers and hospitalization and ICU occupancies. This was accomplished by combining an established compartmental model with a methodology from nonlinear control theory and convolution techniques. No prior assumptions on parameter values have to be made as all parameters can be obtained from data in an online fashion. Thus, users of the method are relieved of the task to find suitable parameter sets of the pandemic model or realistic assumptions for predictions. Predictions up to three weeks into the future for case numbers and occupancies are given for several countries and a statistical assessment quantitatively illustrates the accuracy of the method.
For the online prediction of the aggregated external drivers u(t) simple and robust linear gradient-based functions have been used. Moreover, characteristics and limitations of the course of u(t) have been derived from historic data.
The methodology and its predictive power could be also used to evaluate the effectiveness of measures such as lockdowns or vaccination campaigns, it allows for an early warning system of oncoming waves, and it enables decision makers to evaluate the pandemic situation by incorporating future developments of both state variables and occupancies.
Funding Open access funding provided by TU Wien (TUW).

Data availability
The data used in this paper is taken from the Datenplattform COVID-19 in case of Austria [37], from the publicly available COVID-19 database of data.gouv.fr in case of France [38] and from RIVMdata in case of the Netherlands [39].

Conflict of interest
The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/ by/4.0/.

Estimations of different countries
More prediction results for the age-structured model based on the proposed methodology for different countries are presented in the following: -Predictions of the course of the disease for France are shown in Fig. 11 and for the Netherlands in Fig. 15. -Case-specific admission rates for France are shown in Fig. 12 and for the Netherlands in Fig. 16. -Occupancy predictions for France are shown in Fig. 13 and for the Netherlands in Fig. 17. -Prediction accuracies for France are shown in Fig. 14 and for the Netherlands in Fig. 18.