key: cord-0848366-8gavwbbv
authors: Richard, Q.; Alizon, S.; Choisy, M.; Sofonea, M. T.; Djidjou-Demasse, R.
title: Age-structured non-pharmaceutical interventions for optimal control of COVID-19 epidemic
date: 2020-06-23
journal: nan
DOI: 10.1101/2020.06.23.20138099
sha: fb2cd1855aa73fe4f67378603ad92ea0a13a6ed3
doc_id: 848366
cord_uid: 8gavwbbv

In an epidemic, individuals can widely differ in the way they spread the infection depending on their age or on the number of days they have been infected for. In the absence of pharmaceutical interventions such as a vaccine or treatment, non-pharmaceutical interventions (e.g. social distancing) are essential to mitigate the pandemic. We develop an original approach to identify the optimal age-stratified control strategy to implement as a function of the time since the onset of the epidemic. This is based on a model with a double continuous structure in terms of host age and time since infection. By applying optimal control theory to this model, we identify a solution that minimizes deaths and costs associated with the implementation of the control strategy itself. We also implement this strategy to three countries with contrasted age distributions (Burkina-Faso, France, and Vietnam). Overall, the optimal strategy varies over the course of the epidemic, with a more intense control early on, and depending on host age, with a stronger control for the older population, except in the scenario where the cost associated with the control is low. In the latter scenario, we find strong differences across countries because the control extends to younger population in France and Vietnam 2-3 months after the onset of the epidemic, but not in Burkina Faso. Finally, we show that the optimal control strategy strongly outperforms constant uniform control of the whole population or over its younger fraction. This better understanding of the effect of age-based control interventions opens new perspectives for the field, especially for age-based contact tracing.

Following its emergence in December 2019, COVID-19 has become an international public health emergency [12] . The infection is similar to that caused by influenza virus regarding clinical presentation and transmission mechanism [12] . Contrary to seasonal influenza, COVID-19 has become pandemic by spreading rapidly among completely naive host populations, i.e. with no pre-existing immunity [17, 21, 51, 52] . At the start of the pandemic, no pharmaceutical interventions such as vaccines or treatments were available and, based on earlier epidemics, it will take several months before their deployment. For this reason, developing non-pharmaceutical intervention strategies, such as social distancing, is of great importance to mitigate the pandemic [1] .

Generally, age structure is a key determinant of such acute respiratory diseases, e.g. when it comes to infection severity. For example, children are considered to be responsible for most of the transmission of influenza [8] , but the related hospitalization and mortality burden is largely carried by people of ages over 65 years [38, 50] . While much remains unknown about the COVID-19 epidemics, evidence to date suggests that mortality among people who have been tested positive for the coronavirus is substantially higher at older ages and near zero for young children [40, 51] . Moreover, the infectiousness of an individual has been reported to vary as a function of time since infection [24] , which is known to affect epidemic spread [2, 27, 31] .

Here we propose an epidemiological model for the disease stage-progression [2] structured both by the continuous age of the host population and the continuous age of infection. This formulation differs from the existing literature where only one type of structure is considered at a time [13, 29, 36, 39] , and is particularly suited to investigate an infection such as COVID-19, with strong host and infect age effects. Indeed, in addition to taking into account the host population's age structure, as well as the gradient of disease severity from mild to critical symptoms, the model readily captures the variation in infectiousness as a function of the time since infection. From a theoretical point of view, age-structured models have been proposed to investigate the spread of acute respiratory diseases [4, 18, 34, 47, 48] . However in the literature, very few models consider both structures as continuous variables, see for instance [13, 28] .

In a context of non-pharmaceutical interventions, we adopt a modeling approach based on the optimal control theory to determine the best strategy to implement during a finite time interval. In the context of age-structured models, this approach allows one to determine the optimal strategies of age-specific social distancing taking into consideration the cost of implementing such strategies [3, 5, 6, 14, 25] . Here, more specifically, we look for the intevention that significantly reduces morbidity associated with COVID-19 at a minimal cost. In the same context, mathematical modeling using optimal control theory has been carried out to identify optimal strategies involving non-pharmaceutical interventions to control infectious diseases such as influenza and COVID-19 [15, 30, 35, 41] . However, none of these models take into account the age structure of the host population or the variation of the infectiousness with the time since infection.

In Section 2, we first introduce the mathematical model. The model parameters and outputs are then defined in Section 3. In Section 4, we characterize the optimal control strategy that minimizes the number of deaths as well as the cost due to the implementation of the control strategy itself. Section 5 contains the main body of the results. We first analyse the epidemic spread without any intervention, before comparing the performance of the optimal control in terms of deaths and hospitalizations for 2 different costs of the control measure. Finally, the optimal control is compared to two other strategies using the same amount of resources to control the outbreak. The article ends by a Discussion in Section 6.

2 The age-structured model of COVID-19

At time t ∈ [0, T ], the density of individuals of age a ∈ [0, a max ] that are susceptible to the infection is denoted by S(t, a). These individuals can become infected with a rate called the force of infection and denoted λ (t, a). We assume that a fraction p of these individuals are paucisymptomatic, which means they will develop very mild to no symptoms, and enter group I p . Note that this class can also be interpreted as the fraction of the population that will not isolate themselves during their infection. Other individuals are assumed to develop more symptomatic infections, either severe I s with proportion q(a) depending on the age a, or mild I m with proportion 1 − q(a).

Each of the three infected host populations are structured in time since infection, so that I v (t, a, i), v ∈ {p, s, m}, denotes the density at time t of individuals of age a that have been infected for a duration i ∈ R + . Upon infection, all exposed individuals are assumed to remain non-infectious during an average period i lat . Next, they enter an asymptomatic period during which they are infectious. Only I m and I s develop significant symptoms after an average time since infection i sympt , which can allow them to self-isolate to limit transmission. During their infection, individuals can recover at a rate h v (a, i) (v ∈ {p, m, s}) that depends on the severity of the infection and the time since infection i.

Severely infected individuals of age a may also die from the infection at rate γ(a, i).

The infection life cycle is shown in Figure 1 . The total size of the host population of age a at time t is N(t, a) = S(t, a) + R(t, a) + ∞ 0 (I p (t, a, i) + I m (t, a, i) + I s (t, a, i)) di.

(1)

We use two components to model the infection process. First, we define the transmission probability

for each contact between an infected of age a and a susceptible person, which depends on the time since infection i. Second, we introduce the kernel K(a, a ) that represents the average number of contacts by unit of time between an individual of age a and an individual of age a.

Here, this contact matrix is informed by data from an earlier study conducted in France [7] . The force of infection underwent by susceptible individuals of age a at time t is then given by

S(t,a) hp ( a , i) Figure 1 : The model flow diagram. Susceptible hosts of age a at time t (S(t, a)) are exposed to the virus with a force of infection λ (t, a). A fraction p of exposed individuals, which are infected since time i, will never develop symptoms and enter the group of paucisymptomatic infections (I p (t, a, i)). The rest will develop symptomatic infections, either severe (I s (t, a, i)) with proportion q(a) depending on age a of individuals, or mild (I m (t, a, i)). Exposed individuals remain non-infectious for a duration i lat after infection. Next, they become asymptomatic infectious and only symptomatic infected will develop symptoms at time i sympt after infection. Infected individuals recover at rate h v (a, i). Only severely infected of age a die from the infection at rate γ(a, i). Notations are shown in Table 1 .

Here, c(t, a) is the percentage of reduction of contacts towards people with age a, due to public measures, at time t. The total force of infection at time t in the whole population is computed as 

We assume that only severe infections I s lead to hospitalization and we denote by

the total population hospitalized at time t, where i sympt is the average time to symptoms onset. Each individual of age a dies at a rate µ(a, H(t)) at time t, defined by µ (a, H(t)) = µ nat (a) + µ add (a, H(t)) .

In the latter equation, µ nat denotes the natural mortality rate when hospitals are not saturated. Further, we assume that this rate increases significantly as soon as the number of severe cases exceeds the healthcare capacity H sat and µ add is such additional death rate due to hospital saturation (see Section 3.2). We apply the same reasoning by assuming that the disease-related mortality can increase because of hospital saturation. Therefore, severely infected individuals of age a infected since time i die at time t at rate γ(a, i, H(t)) defined by γ(a, i, H(t)) = (γ dir (a) + γ indir (a, H(t))) 1 [i sympt ,i s max ] (i).

Here, γ dir and γ indir are mortality rates directly and indirectly due to the COVID-19 respectively (see Section 3.2). The disease-related mortality occurs after the emergence of symptoms and before the mean final time of infection for severe cases, i.e. for i ∈ [i sympt , i s max ]. Finally, infected individuals of age a infected since time i recover at rates h s (a, i), h m (a, i) and h p (a, i) for severe, mild and paucisymptomatic infections respectively.

The boundary conditions (3) are coupled with the following equations:

for any (t, a, i) ∈ (0, T ] × [0, a max ] × R + , with initial conditions (at t = 0):

The initial conditions of infected populations are detailed in Section 3.3. Using (3) and an integration over i of (5), one may observe that the total population N defined by (1) is strictly decreasing since it satisfies the following inequality:

This is due to the fact that population aging and births are neglected in this model since we consider a time horizon of only one year.

In this section we briefly describe some useful epidemiological outputs, the shape of age dependent parameters considered for the simulations of model (3)-(5), and the initial conditions. All state variables and other parameters are summarized in Table 1 .

In addition to the total number of hospitalized cases H(t) at time t defined by (4), we define additional epidemiological outputs such as the number of non-hospitalized cases (N H (t)), the cumulative number of deaths due to COVID-19 directly (D cum dir (t)) and indirectly (D cum indir (t)) respectively by

5 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted June 23, 2020. 

where D dir (t) and D indir (t) are the number of deaths at time t respectively defined by

Note that non-hospitalized cases N H defined by (6) are composed of the paucisymptomatic, the mildly infected, and the severely infected but not yet hospitalized populations. We can also note that every output aforementioned implicitly depends on parameter c which we will omit when no confusion is possible. However, in order to compare different public health measures we will explicitly write this dependence. The relative performance between two strategies c 1 and c 2 , denoted by ∆(c 1 , c 2 ), is estimated by assessing the cumulative number of deaths in the whole population during the T days of control period with the strategy c 1 relatively to deaths with the strategy c 2 . Formally we have

Hence, a relative performance ∆(c 1 , c 2 ) = 0.1 implies that the strategy c 1 reduces the number of deaths by 10% relatively to c 2 .

We assume mortality rates indirectly due to the COVID-19 to grow as the number of hospitalisations H exceeds a healthcare capacity threshold H sat . The natural mortality rate increases by µ add (a, H) in the whole population and by γ indir (a, H) for severely infected individuals of age a. These rates are respectively defined by logistic functions that are arbitrarily chosen as:

This choice of functional parameters implies that

so that those additional mortalities are negligible when hospitals are not saturated (Figure 2 b,c). In case of saturation, the following estimates hold:

for each a ∈ [0, a max ], meaning that the natural mortality rate is only increased by 1%, while the disease-induced mortality rate γ is doubled. Indeed, according to [23], 50% of patients in critical care 6 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. will die in case of no saturation of hospitals. Here, we then make the assumption that this percentage will grow to 100% in case of over-saturation of hospitals.

The infectiousness of an individual aged a, which is infected since time i, is given by β v (a, i) (v ∈ {s, m, p}). Based on estimates described in [24], we assume that β v does not depends on age a, i.e.,

This assumption is discussed later in Section 6. Next, we set

Here, as explained below, α is a scaling parameter obtained from the value of the basic reproduction number R 0 . Parameter β is assumed to be identical to that reported in [24] and to follow a Weibull distribution β ∼ W (3, 5.65). Parameters ξ v (i) are factors capturing the reduction of the transmission probability. For paucisymptomatic individuals, these are assumed to be constant (ξ p (i) = ξ p ), while the reduction factor in more symptomatic infections (severe and mild) is assume to vary after symptom onset to capture admission in a healthcare facility or self-isolation at home. More precisely, we assume that

These two functions are chosen arbitrarily by assuming that individuals do not isolate before symptoms onset (i ≤ i sympt ), and that isolation is stronger when symptoms are more ( Figure 2 (a)). We therefore assume that the transmission probability β is divided by 10 (respectively 2) every day after the average time of symptoms onset for individuals severely (resp. mildly) infected. Finally, we assume that recover rates h v (a, i), v ∈ {s, m, p}, of infected individuals of age a infected since time i are independent of the age a and take the following form:

That is, one can recover from severe (resp. mild and paucisymptomatic) infections only after a time since infection i s max (resp. i m max ) corresponding to the mean duration of infection.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint

According to the French public health agency [16] , there were 130 confirmed cases of COVID-19 in France on March 1st, 2020, which we will consider as t = 0 in our model. Since tests in France were initially performed based on severe symptoms, we assume that all those cases are severe infections. Thus, we set i s max i sympt a max 0 I s,0 (a, i)da di = 130 as the initial severely infected individuals, which is assumed to be uniformly distributed with respect to the time since infection i on the interval [0, i s max ]. Using estimates from [10, 16] on the age distribution of hospitalised people, we derive an estimation of I s,0 (a, i) for each (a, i) ∈ [0, a max ] × R + . Next, following the life cycle (Figure 1) , we obtain an estimation of the total initial infected population by

. From there, we deduce the initial mildly and paucisymptomatic infected populations respectively by

The initial susceptible population size S 0 comes from the French National Institute of Statistics and Economic Studies [20].

In this section, following well established methodology in optimal control theory [3, 6, 14, 22 , 25], we search for the optimal control effort function c * that minimizes the objective functional J :

where D cum dir , D cum indir are cumulative number of deaths defined by (7) and B(a) is the cost associated with the implementation of such control c for the age class a. Our aim is to find the function c * satisfying

wherein the set U is defined by

with c max ≤ 1 a positive constant. That is to say, the function c * will minimize the cumulative number of deaths during T days, as long as the cost of the control strategy is not too large. Let (S, I s , I m , I p , R) be a given solution of (3)-(5) then let λ and H be respectively defined by (2) 9 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint and (4). After some computations (Appendix B), we find that the adjoint system of (5) reads as

with final conditions z S (T, a) = z R (T, a) = 0, z u (T, a, i) = 0 and lim i→∞ z u (t, a, i) = 0, for any u ∈ {I s , I m , I p } and (a, i) ∈ [0, a max ] × R + , while ζ k (k ∈ {1, 2, 3}) satisfy the system:

Finally, the Hamiltonian H of (11) is given by (B.1). Then, solving ∂ H ∂ c = 0, it comes that c * (t, a) = max(0, min(ĉ(t, a), 1)),

for 

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. .

where B * ∈ R + is a variable parameter characterizing the relative cost in implementing the strategy.

The state system (3)-(5) and the adjoint system (12)- (13) together with the control characterization (14) form the optimality system to be solved numerically. Since the state equations have initial conditions and the adjoint equations have final time conditions, we cannot solve the optimality system directly by only sweeping forward in time. Thus, an iterative algorithm, forward-backward sweep method, is used [32] . In other words, finding c * numerically, involves first solving the state variables (3)-(5) forward in time, then solving the adjoint variables (12)-(13) backward in time, and then plugging the solutions for the relevant state and adjoint variables into (14) , subject to bounds on the control function. Finally, the proof of the existence of such control is standard and is mostly based on the Ekeland's variational principle [19] . Therefore, existence of the optimal control to the above problem is assumed and we refer to [14] for more details.

An explicit expression of the R 0 of model (3)-(5) is difficult to obtain in general. We show in Appendix A that it is possible to write R 0 = α × r(U), where α is the scaling parameter introduced in Section 3.2, and r(U) is the spectral radius of the next generation operator U defined on L 1 (0, a max ) by

where S 0 is the initial susceptible population, K is the contact matrix and ω(a, i) is the infectiousness of individuals of age a infected since time i (Appendix A). We set R 0 = 3.3 [45, 49] and it follows that

Using a numerical approach, we find r(U) ≈ 12.074, whence α ≈ 0.575.

Numerical simulations are based on the reference values of the model parameters defined previously and summarized in Table 1 , with R 0 = 3.3.

We first use the model to describe the outbreak of the epidemics without any public health measure (i.e. c ≡ 0). The peak of the epidemics is reached approximately at day t = 54 for hospitalised people, and day t = 48 for non-hospitalised (Figure 3a) . The delay between the two peaks is due to the latency time i sympt for symptoms onset ( Table 1 ). The healthcare capacity is quickly exceeded (about twenty days) and the number of deaths increases sharply from then on. At the end of the simulation (t = 150 days), the total number of infections (severe, mild and paucisymptomatic) is around 90.1% and varies with the age class considered (Figure 3 b) . Further, in each age class the proportion of 11 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint infected individuals is larger than the theoretical herd immunity threshold given by 1 − 1/R 0 ≈ 69.7% (Figure 3b ). While people older than 70 are the less affected (in proportion), they however represent the age class with the highest cumulative number of deaths (Figure 3c ). On the contrary, most of the infections that occur in the young population do not require hospitalisation (Figure 3d ,e).

In this section, we investigate the interaction between the optimal intervention and the age structure of the population. We illustrate the optimal intervention strategy and their performance in terms of cumulative number of deaths for three costs of the control (relatively low B * = 10 2 , intermediate B * = 10 3 and high B * = 10 4 ). Overall, the optimal control particularly targets the older populations compared to the younger ones ( Figure 4 ). If the cost B * is relatively high, the optimal control is almost restricted to individuals above 55, with a significant reduction in deaths (Figure 4d ,e). This strict lock down of older individuals lasts approximately 100 days for B * = 10 3 and 10 4 . The relative performance of the optimal control c * compared to a doing nothing scenario (∆(c * , 0)) is at least 92% (resp. 82%) when the cost is B * = 10 3 (resp. 10 4 ). With a low cost of the control measure (B * = 10 2 ), the optimal control significantly extends to younger populations (Figure 4 a) , with a maximum reached near the 4th month of the epidemics and a steady decrease until the end of the control period. At first, the control is over people above 35 but that after 2 or 3 months the control begin to extend to people less than 35. The resulting reduction in the number of deaths is very pronounced with a relative performance ∆(c * , 0) of at least 99%.

We investigated how the optimal strategy compares to other control strategies that use the same amount of 'resources' (that is the same cumulative cost). Assuming an relatively high cost B * = 10 3 , we investigated a uniform control strategy (denoted by c u ) either over the younger fraction of the population (Figure 5a ) or over the whole population with level c max = 0.95 (Figure 5b ). The effect of these strategy lasts about 55 days during which the epidemic is suppressed. However, once the control resources are exhausted, the epidemic reemerges ( Figure 5c ) and, in the end, the cumulative mortality over the time period of interest is comparable to that without any control measure ( Figure 5d ) and relative performance ∆(c * , c u ), of the optimal control c * relatively to the uniform control c u , is approximately 92%. With the (longer) uniform control over the younger fraction of the population, the first epidemic peak is slightly delayed, but a second peak appears a few months later (Figure 5c ). With this strategy, the cumulative mortality comparable to the one without any control measure (Figure 5d ).

The optimal control is a continuous function and is then quite difficult to enforce in practice. However, we can derive step functions leading to practical implementations of the optimal control. For instance, the population is subdivided into 10-year amplitude classes and the control is updated every 3-weeks by keeping a constant amount of control during each 3-weeks period for each age-class.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. 13 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint ( Figure S1 ). Importantly, the constant defining the control intensity for each period is captured from the knowledge of the continuous optimal control strategy. The effect of such strategy is overall similar to the optimal control ( Figure S1 ) with a relative performance of 91% compared to a doing nothing scenario. (e) Cumulative deaths per age at final time T = 365 days. The relative performance ∆(c * , 0) of the optimal control c * compared to a doing nothing scenario is at least 99% (resp. 92%, 82%) with B * = 10 2 (resp. 10 3 , 10 4 ). Here p = 0.5.

Non-pharmaceutical public health interventions can be implemented either to mitigate the COVID-19 epidemic wave, and rely on natural immunisation, or to suppress the wave long enough to develop and implement a vaccine or a treatment. Here, we explicitly factor in the age heterogeneity of the host population in the identification of the optimal allocation of the control efforts. We use optimal control theory to characterize an optimal strategy that significantly reduces the number of deaths, while being sustainable at the population level. Our formulation assume a quadratic cost for the control effort. We find that, with this strategy, the intensity of the control is always relatively high on the older fraction of the population during at least a hundred days, before decreasing 14 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. 

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint more or less rapidly depending on the cost associated to the control. The control over the younger fraction of the population is weak and only occurs when the cost associated with the optimal control is relatively low and, even then, the level control only increases 2 or 3 months after that on the older fraction of the population. This late control over the younger part of the population actually mimics the results [15] where the control didn't peak right away. Intuitively, if control strategies come at a high cost for the population, it is best to focus on the age classes that are the most at risk. Conversely, if the control measures are more acceptable to the population, the optimal strategy is to aim wide in order to completely suppress the epidemic wave.

Information on the natural history of paucisymptomatic infections of COVID-19 remains relatively little-known [9, 44] . It is estimated that a proportion p of infected individuals will remain asymptomatic throughout the course of infection. However, this proportion remains largely unspecified in the literature [9, 44] . We explored effects of the proportion p on the optimal control strategy c * . Overall, the proportion of paucisymptomatic infections have marginal effect of the optimal control strategy ( Figure S2 ). The optimal control remains strong over the older population from the beginning of the epidemic, before progressively alleviated. The control over the younger population is weaker and occurs only if control cost themselves is low. But, the level of control over its younger fraction increases when the proportion of paucisymptomatic infections decreases.

Given the leverage represented by school and university closure, we investigated the effect of control measures over individuals aged under 25. Our results show that NPIs targeting the younger fraction of the population are not very efficient in reducing cumulative mortality, unless they can be implemented strongly and for a relatively long period. Indeed, a control only over the younger population barely reduces the number of deaths by 3% compared to a doing nothing scenario (Figure 5a,c) . Note that this result could depend on the contact matrix between ages, for which there is little date in France (here we used the one found by [7] ). Furthermore, transmission probability could vary with age, as discussed below.

The model proposed here is an extension of the classical models based on ordinary differential equations that tackled the issue of the optimal control of COVID-19 outbreak [15, 30, 35, 41] . Here, the whole population is structured by age (a) and additionally by the time since infection (i) for infectious individuals, which echoes the model developed in [49] using a discrete-time formulation of the infection. With our continuous structure, we show that the number of new cases I N (t, a) at time t in individuals of age a is given by the renewal equation

where K is the contact matrix and ω(a, i) is the infectiousness of individuals aged a which are infected since time i (Appendix A). For parameterisation purpouses, we assume that ω(a, i) is the product between the proportion of individuals of age a in the whole population and the infectiousness β (i) of individuals infected since time i. This is potentially a strong limitation since infectiousness β could depend on the age a thereby creating an additional heterogeneity in addition to that since the time 16 . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . since infection i. This issue can be particularly important since some studies suggest a low risk of transmission in the young population (e.g. [11] ).

Another potential limiation is the lack of gender structure in the model formulation. Given the observed male biased in mortality during the COVID-19 pandemic, it has been suggested that males are more at risk of developing severe infections [46] . This heterogeneiry could readily be introduced in the model.

Contacts networks have an important role in transmission dynamic models. Epidemic models that determine which interventions can successfully prevent an outbreak may benefit from accounting for social structure and mixing patterns. Contacts are highly assortative with age across a given country, but regional differences in the age-specific contacts is noticeable [42] . The current model could be modified to explore epidemiological dynamics in a spatially structured population with nonhomogeneous mixing, e.g. by using a meta-population model [37] .

Another potential extension of the model would be to allow for the isolation of symptomatic cases and their contacts, following the method developed in [26] and applied recently to digital contact tracing [24] . Indeed, these measures strongly depend on the relative timing of infectiousness and the appearance of symptoms, and the formulation of the presented model seems suitable for that. However, this also raises technical challenges due to the double continuous structure. However, being able to identify age classes to follow in priority with contact tracing could be an asset in controling epidemic spread.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint . CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. [25] Fister, K. R. and S. Lenhart. 2004 

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. Figure S1 : Practicability of the age-structured optimal control. (a)-(b)

Step optimal controls with a 3-weeks update over the older and younger populations. The corresponding optimal is given by Figure 4b . (c)-(d) Cumulative deaths per age at final time T = 365 days.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. 

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint

Here we compute the basic reproduction number R 0 of the model (3)- (5) . First let us set for i ≥ 0 and a ∈ [0, a max ] the following functions

that describe the survival probability of infected individuals (in the respective compartment), with age a, from their infection until the time since infection i, in case of no hospitalisation (i.e. H ≡ 0). We get the following Volterra formulation of the linearized system of (3)- (5): 

and f (t, a) is the density of new infections produced by the initial population. Therefore, the basic reproduction number R 0 is the spectral radius, denoted by r(U), of the next generation operator U defined on L 1 + (0, a max ) by

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint As explained in Section 3.2, it is estimated in [24] that each average infectiousness β k (k ∈ {s, m, p}) takes the form of a Weibull distribution W (3, 5.65) so that the mean and median are equal to 5.0 days while the standard deviation is 1.9 days. Based on this estimation, we assume that β k (a, i) = αβ (i)ξ k (i) where β ∼ W (3, 5.65) and α is a positive parameter to be determined. Consequently, it follows that α is given by

where U is the operator defined by

We see that U can be rewritten as Now, in order to compute the spectral radius r U , we first make the following assumptions:

Assumption A. 1 We suppose that: a) functions S 0 , K, Ω are bounded and positive almost everywhere; b) the function K is symmetric.

We can note that the Assumption A.1 is satisfied when using the parameters stated in Table 1 . Now, let S be the positive self-adjoint operator defined by S : L 2 (0, a max ) v −→ S 0 (·)Ω(·) Proof. The compactness of both integral operators follows from the fact that a max < ∞ by assumption (see Table 1 ), hence their spectra are punctual. Now we prove that σ (U) = σ (S). Let ν ∈ σ (U) be an eigenvalue of U and φ ∈ L 1 (0, a max ) be the associated eigenvector, i.e.

Uφ (a) = S 0 (a) 

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint so that φ ∈ L ∞ (0, a max ) ⊂ L 2 (0, a max ). Defining the function i.e. ν ∈ σ (S) is an eigenvalue of S associated to the eigenvector ψ, so that σ (U) ⊂ σ (S). For the reverse inclusion, let ν ∈ σ (S) and ψ ∈ L 2 (0, a max ) ⊂ L 1 (0, a max ) be the associated eigenvector for S. It follows that the function

is an eigenvector of U related to the eigenvalue ν ∈ σ (U), whence σ (U) = σ (S). In particular, both spectral radius are equal. Finally, the Rayleigh formula is classical for positive and symmetric operators.

Remark A.3 Numerically, to compute r(U), we can similarly show that it is given by the spectral radius of the following operator: which can be easily computed since the age a is numerically divided into 20 classes, so that the term inside the integral of the latter equation is a 20 × 20 matrix. Finally, we obtain α from (A.5).

In order to deal with the necessary optimality conditions, we use some results in [22] . Next, we detail the computations of the adjoint system (12)- (13) . To this end, we first define the functions g R (i, y 2 (t, a, i))di.

. CC-BY-NC 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted June 23, 2020. . https://doi.org/10.1101/2020.06.23.20138099 doi: medRxiv preprint

The model (5) thus rewrites as      ∂ t y 1 (t, a) = F 1 (a, Q(t, a), c(t, a), y 1 (t, a)), (∂ t + ∂ i ) y 2 (t, a, i) = F 2 (a, i, Q(t, a), c(t, a), y 2 (t, a, i)), y 2 (t, a, 0) = Φ(a, c(t, a), E(t, a)), with F 1 (a, Q(t, a) , c(t, a), y 1 (t, a)) = −µ(a, H(t))S(t, a) − (1 − c(t, a))E(t, a) −µ(a, H(t))R(t, a) + b(t, a) ,

Special Report: The Simulations Driving the World's Response to COVID-19

Infectious Diseases of Humans. Dynamics and Control

Analysis and Control of Age-Dependent Population Dynamics

Transmission Dynamics of Acute Respiratory Diseases in a Population Structured by Age

Optimal Intervention Strategies of Staged Progression HIV Infections through an Age-Structured Model with Probabilities of ART Drop Out

Optimal Control of Population Dynamics

The French Connection: The First Large Population-Based Contact Survey in France Relevant for the Spread of Infectious Diseases

Identifying Pediatric Age Groups for Influenza Vaccination Using a Real-Time Regional Surveillance System

The Role of Asymptomatic SARS-CoV-2 Infections: Rapid Living Systematic Review and Meta-Analysis

Severe Outcomes Among Patients with Coronavirus Disease 2019 (COVID-19) -United States

Assessment of Spread of SARS-CoV-2 by RT-PCR and Concomitant Serology in Children in a Region Heavily Affected by COVID-19 Pandemic

Coronavirus Disease (COVID-19) Situation Reports (n.d.)

Proportionate Mixing Models for Age-Dependent Infection Transmission

Optimal Control for an Age-Structured Model for the Transmission of Hepatitis B

Optimal COVID-19 Epidemic Control until Vaccine Deployment

Données hospitalières relatives à l'épidémie de COVID-19 -data.gouv.fr

2020. Beyond Just "Flattening the Curve": Optimal Control of Epidemics with Purely Non-Pharmaceutical Interventions

A Contribution to the Mathematical Theory of Epidemics

Optimal Control Applied to Biological Models

Leung and Z. Feng. 2020. Early Transmission Dynamics in Wuhan, China, of Novel Coronavirus-Infected Pneumonia

Deep Reinforcement Learning for Large-Scale Epidemic Control

An Optimal Control Theory Approach to Non-Pharmaceutical Interventions

Predicting the Number of Reported and Unreported Cases for the COVID-19 Epidemic in South Korea

Spatial Heterogeneity and the Design of Immunization Programs

New Estimates of Influenza-Related Pneumonia and Influenza Hospitalizations among the Elderly

Global Stability for an SEI Epidemiological Model with Continuous Age-Structure in the Exposed and Infectious Classes

Case-Fatality Rate and Characteristics of Patients Dying in Relation to COVID-19 in Italy

Projecting Social Contact Matrices in 152 Countries Using Contact Surveys and Demographic Data

Pyramide Des Âges | Insee

Natural History of Asymptomatic SARS-CoV-2 Infection

Considering How Biological Sex Impacts Immune Responses and COVID-19 Outcomes

Optimal Strategies of Social Distancing and Vaccination against Seasonal Influenza

Age-Structured Impact of Social Distancing on the COVID-19 Epidemic in India

Epidemiological Monitoring and Control Perspectives: Application of a Parsimonious Modelling Framework to the COVID-19 Dynamics in France

Influenza-Associated Hospitalizations in the United States

Estimating Clinical Severity of COVID-19 from the Transmission Dynamics in Wuhan, China