key: cord-0816106-vz7do8us authors: Campos, Eduardo Lima; Cysne, Rubens Penha; Madureira, Alexandre L.; Mendes, Gélcio L.Q. title: Multi-generational SIR modeling: Determination of parameters, epidemiological forecasting and age-dependent vaccination policies date: 2021-06-10 journal: Infect Dis Model DOI: 10.1016/j.idm.2021.05.003 sha: d08291e6592eacb0749abdd584fe7f9f93a93f9f doc_id: 816106 cord_uid: vz7do8us We use an age-dependent SIR system of equations to model the evolution of the COVID-19. Parameters that measure the amount of interaction in different locations (home, work, school, other) are approximated from in-sample data using a random optimization scheme, and indicate changes in social distancing along the course of the pandemic. That allows the estimation of the time evolution of classical and age-dependent reproduction numbers. With those parameters we predict the disease dynamics, and compare our results with out-of-sample data from the City of Rio de Janeiro. Finally, we provide a numerical investigation regarding age-based vaccination policies, shedding some light on whether is preferable to vaccinate those at most risk (the elderly) or those who spread the disease the most (the youngest). There is no clear upshot, as the results depend on the age of those immunized, contagious parameters, vaccination schedules and efficiency. Several months after the onset of the COVID-19 pandemic it became clear that good numerical modeling is useful to understand the dynamics of the disease. In particular, simple compartmental epidemiological models are convenient, depending on few parameters and still capturing with reasonable precision essential aspects of the dynamics of infectious diseases (Brauer, 2017; Hethcote, 2000; Wang et al., 2021) . The Susceptible-Infected-Recovered (SIR) model in its simplest form depends on two parameters, the average number of sufficient contacts (those sufficient for transmission) and the mean waiting period during which patients are contagious. It relies on the assumption that the population is homogeneous both in terms of behavior and biological susceptibility. It also assumes, in its simplest instance, that demographic aspects as age structure, birth and death rates do not matter for the dynamics of the disease. Some of the above assumptions appear to hold well in the case of the COVID-19 pandemic. For instance, demographic changes appear to be of less importance due to difference in time scales, compared with the evolution of the pandemic. However, dropping the homogeneity hypothesis allows a more realistic model, as age and location of contacts are important factors in understanding the dynamics of the disease. For instance, the disease's severity and lethality is age-dependent, among other factors (Bhopal & Bhopal, 2020; Jin et al., 2020; Verity et al., 2020; Wu & McGoogan, 2020) . Mathematical epidemiology has a somewhat long history, with (Hamer, 1906; M'Kendrick, 1925; Kermack, McKendrick, & Walker, 1927) being among the earliest contributions. See (Brauer, 2017; Hethcote, 2000) for a thorough review of mathematical modeling of infectious diseases, and (Albani et al., 2021a Iannelli & Milner, 2017; Klepac & Caswell, 2011; Klepac et al., 2009; bib_Li_et_al_2020Li, Yang, & Martcheva, 2020 Prem et al., 2020; Sun, 2010; Towers & Feng, 2012) for multi-group models. Some other models consider spatial aspects of infectious diseases, as (Bertaglia & Pareschi, 2020; Besse & Faye, 2020; Lang, De Sterck, Kaiser, & Miller, 2018; Paeng & Lee, 2017; Peixoto, Marcondes, Peixoto, & Oliva, 2020; Tak acs & Hadjimichael, 2019; Viguerie et al., 2020) . Theoretical aspects of models are covered in (Kuniya & Wang, 2018; Wen, Ji, & Li, 2018; Wu & Zou, 2016) . For the control related models for COVID-19 see for instance (Perkins & España, 2020) and references therein. After the onset of the COVID-19 pandemic, a massive amount of articles related to the topic materialized. Several articles considered compartmental models to forecast different scenarios (Albani et al., 2021a Amaku et al., 2021a , bib_Amaku_et_al_2021a , bib_Amaku_et_al_ 2021b Atkeson, 2020; Coelho et al., 2020; Kimathi, Mwalili, Ojiambo, & Gathungu, 2021; Komatsu & Menezes-Filho, 2020; Prem et al., 2020; Viguerie et al., 2020; Walker et al., 2020) . Socioeconomic costs of imposing social distancing are combined with compartmental models in (Acemoglu, Chernozhukov, Werning, & Whinston, 2020; Alvarez, Argente, & Lippi, 2020; Atkeson, 2020) , and the interplay between economics and individual decisions are considered in (Borelli & G oes, 2020; Brotherhood, Kircher, Santos, & Tertilt, 2020; Eichenbaum, Rebelo, & Trabandt, 2020) . See (Kremer, 1996; d'Onofrio & Manfredi, 2020, pp. 185e203; Funk et al., 2015; Manfredi and d'Onofrio, 2013; Perrings et al., 2014; Soofi, Najafi, & Karami-Matin, 2020 ) for a description of behavioral-and economic-epidemiology ideas. Also, macro-economic aspects of the pandemic appears, e.g., in (Guerrieri, Lorenzoni, Straub, & Werning, 2020) . Population heterogeneity and its influence on herd immunity is considered in (Aguas et al., 2020; 2020Britton, Ball, & Trapman, ; Cui, Zhang, & Feng, 2019; Gomes et al., 2020) . Also, different aspects of vaccination are considered in the literature (Albani et al., 2021b; Amaku et al., 2021a; d'Onofrio & Manfredi, 2020, pp. 185e203; Meehan et al., 2020; Chen & Toxvaerd, 2014; Gonz alez & Villena, 2020; Jia et al., 2020; May & Anderson, 1984) . Finally, alternative approaches, as agent-based, statistical models and a combination of these techniques with compartmental models (specially to gather data), are also commonly employed (Amaral, Casaca, Oishi, & Cuminato, 2021; Bendtsen Cano et al., 2020; Calvetti, Hoover, Rose, & Somersalo, 2020; Donnat & Holmes, 2004; Ferguson, Laydon, & Nedjati-Gilani, 2020; Roda, Varughese, Han, & Li, 2020) . In the present paper, we consider a modification of the age-structured SIR model (Britton et al., 2020; Towers & Feng, 2012) , in particular allowing for age-dependent probabilities of patients being sub-clinical (asymptomatic or symptomatic patients that does not require hospitalization) or clinical (Albani et al., 2021b; Prem et al., 2020) . We use contact matrices (Prem, Cook, & Jit, 2017) which differentiates not only age but also location (home, work, school, other); see also (Kimathi et al., 2021) . Random searches approximate some key parameters that are related to social distancing in different locations, computing their values for some Brazilian regions. Following the approximation of those parameters, we are able to forecast long term best/worst case scenarios. We gauge the efficiency of the scheme comparing our predictions with available data using in-sample/out-of-sample testing. We next use the model to compare age-dependent vaccination policies (Meehan et al., 2020) . In a simple setting, we compare vaccinating the young or the older population, since the model allows the immunization of specific age groups. The outcome depends on the values of certain parameters, including how aggressive is the transmission and how late after the start of the pandemic the vaccinations take place. In general, for a short term period, vaccination prioritizing the elderly population is the best option with respect to the total number of deaths. But that might not be the case in the long run if extra vaccine supplies are not provided. Disclaimer: It is not our intention to propose real-life epidemiological policies or to provide predictions upon which policies could be developed. Our intention is solely to investigate possible patterns of the disease evolution and raise possible alternatives in terms of vaccination policies. In the spirit of compartmental modeling (Hethcote, 2000; Wang et al., 2021) , we divide the whole population among susceptible (S), infected (I), removed (R). Due to the short time scale of the disease, we can assume that demographic changes are not relevant. The number of deaths are estimated a posteriori, from the removed R. The above quantities are 16-E.L. Campos, R.P. Cysne, A.L. Madureira et al. Infectious Disease Modelling 6 (2021) 751e765 dimensional vectors, stratified by age, running from 0 to 4 (first components), 5e9 (second components), etc, up to 75-above (16th components). There are P i individuals for each age group, and the total population P ¼ P 16 i¼1 P i is constant. The equations governing the dynamics of the disease is as follows, 1 for i ¼ 1 …, 16 and t2N: (1) plus initial conditions. The model is "conservative" in the sense that S i þ I i þ R i ¼ P i is constant for all ages i. The term b i S i is the rate of new infections of the i-population, and we remark that b depends on the number of patients of all ages, making the model coupled and nonlinear. These patients are removed from their condition at rate gI i (they either recover or die). We further divide the infectious into two groups, those that are asymptomatic or sub-clinical I sc i ¼ ð1 Àr i ÞI i in the sense that they do not require medical care. The clinical group I c i ¼ r i I i corresponds to those that get ill, requiring medical attention. The age-dependent fractions of those who become sub-clinical or clinical are r i 2 [0, 1]. To fully describe the model, we ought to define b and the other parameters. (i) the components of the incidence function b are given by the formula where C e is an effective contact matrix that measures the number of contacts between age groups. It depends on the amount and choices of lock-down imposed, as we describe further ahead. The age-independent parameters a sc and a c correspond to the fraction of those that are infective, having the potential to infect others due to social behavior or biology. For convenience, we also write the vector b(I) ¼ C e DI, where D is a diagonal matrix defined by D jj ¼ a sc ð1 À r j Þ þ a c r j P j : (ii) The time dependent effective contact matrix C e is defined by where the exogenous data C home , C work , C school , C other are 16 Â 16 matrices, where C [ ij indicates the average number of different people from the age group j that someone from the age group i contacted per day at [ 2{home, work, school, other}, as compiled by (Mossong et al., 2008; Prem et al., 2017) . The parameters b [ control the fraction of contacts that are "adequate", in the sense that they lead to infections. (iii) Each component of r gives the risk that an infectious will become a sub-clinical or a clinical case. Since symptoms are more aggressive for older patients, r i grows with i. (iv) g ¼ 1 À e À1=dI is the daily probability that someone stays infected, where d I is the average duration of infection (Prem et al., 2017) . For simplicity, we assume that the period in which one is infective is independent of age and disease severity. These hypotheses are commonly, but not always, assumed in the literature (Acemoglu et al., 2020; Albani et al., 2021a Britton et al., 2020; Prem et al., 2020) ; see (Levin et al., 2020; Walsh et al., 2020; Zhou, She, Wang, & Ma, 2020) for related discussions based on medical data. (v) The contagious parameters a sc and a c correspond to fractions of sub-clinical and clinical infected patients that still spread the virus. The basic reproduction number R 0 is understood as the secondary cases produced by an infective patient introduced in a totally susceptible population (1002Barril, Calsina, Cuadrado, & Ripoll, ; Delamater, Street, Leslie, Yang, & Jacobsen, 2019; Diekmann, Heesterbeek, & Metz, 1990; Hethcote, 2000) . Hence, suppose that the whole population with distribution P is susceptible (S ¼ P), and we represent the infective patients byÎ, a vector belonging to the canonical basis {e 1 , …, e 16 }. The number of patients in the first day belonging to the i-th group is where d ik denotes the Kronecker delta, which equals one if i ¼ j and zero otherwise, and we define P as the diagonal matrix Then the vector containing the contaminated individuals in the first day is AÎ, and we can bound its Euclidean norm as where l max is the largest eigenvalue of A. Note that only the fraction ð1 À gÞ dÀ1Î of the infective outsider remains infected on the d-th day. This fraction infect less than (1 À g) dÀ1 l max , as in (3), and adding up all days, we gather that the total number of infective agents is bounded by the basic reproduction number Reasoning in a similar fashion, we define the replacement number where l max (t) is defined as the maximum eigenvalue of the matrix SC e D and S is the diagonal matrix with the diagonal given by S j for j ¼ 1, …, 16. The above definitions of R 0 and R t rely on the Euclidean and spectral norms. However, if one is interested in considering age-dependent quarantine or vaccination policies, it might be interesting to consider a slightly different approach. An infective patient of the age group [ 2{1, …, 16} would contaminate, on a single day, As before, adding up the days and considering that at the jth day there will be a fraction (1 À g) dÀ1 e [ of the infective outsider, we gather that the total number of infective agents is the age-dependent basic reproduction number Similarly, we define the corresponding age-dependent replacement number The above number R t;[ indicates how many people would be contaminated by an infective patient from the age group [ introduced in the community at time t. Mortality is not part of the dynamics, 2 resulting from the number of infectious. We postulate that only those clinically infected die, at an age-dependent rate. So, among the "recovered" patients, die, where the constant vector w d is the age-dependent lethality weights of the disease (Levin et al., 2020) , and the lethality strength m is a scalar that captures how lethality varies with time. Note that the lethality proportion between different ages is constant. 2 We assume that the recovered patients do not get reinfected. The transmission parameters a sc and a c are such that the transmission rate of symptomatic is higher than those asymptomatic (Byambasuren et al., 2020; Prem et al., 2020) , and we stipulate that the proportion of symptomatic patients follows the same pattern of the hospitalized patients in . However, the proportion of asymptomatic patients is far from being settled (Byambasuren et al., 2020; Heneghan, Brassey, & Jefferson, 2020) , with estimates ranging from 5% to 80%. Also, we assume that profile of the lethality rates follow Fig. 1 , based on (Brazeau et al., 2020; Verity et al., 2020) . The contact matrices C home , C work , C school , C other are from (Prem et al., 2017) , and population data, as seen in Fig. 10 , is from the Instituto Brasileiro de Geografia e Estatística (IBGE). The average number of days d I that a patient is infective is set to 12. We use the available data (number of new cases and deaths) from each location to estimate the lethality strength m(t), , which are necessary to compute the effective matrix C e . These functions are transient since not only individuals change behavior with time (Brotherhood et al., 2020) , but also lethality changes with medical conditions and virus mutations. For a simple SIR model, the available data uniquely determines the parameters (Bakhta, Boiveau, Maday, & Mula, 2020) . In general, that is not the case for our model, for instance if one of the four contact matrices vanish or if C work ¼ C other . We reduce the possibility of parameter indeterminations by assuming that b w ¼ b o , i.e., the behavior at work is equal to the behavior at other locations. We also make then some further simplifying hypotheses. First, we assume that b h is constant, as there is no reason to expect that the amount of contacts at home varies significantly with time. The second assumption is that b s (t) is identically zero, since schools were closed during the period of data gathering. Consider that data is available for the period f1; …; Tg4N, and let and N data I ðjÞ and N data D ðjÞ be the number of new cases and new deaths at the j-th day. For simplicity, assume that T ¼ Ndt where dt ¼ 10, and consider the partition T dt of f1; …; Tg Let P 0 ðT dt Þ be the space of piecewise constant functions with respect to the above partition. To determine b w we pose the problem of finding (b h , b w ) that solves and k ,k R T is the Euclidean norm in R T . Above, N sir I ðb h ;b w ÞðtÞ is the number of new infected clinical patients that the SIR model (1) yields at the t-th day if one replaces b h and b w byb h andb w in (2). The above problem has finite dimension N þ 1 since any function b2P 0 ðT dt Þ can be written as b w ðtÞ ¼ for some unknown coefficients b 1 ;…;b N 2R. Here, c Ij ðtÞ is the characteristic function of I j , which equals one if t 2 I j , and zero otherwise. Note that the coefficients uniquely define a function in P 0 ðT dt Þ and vice-versa. To solve (10) we employ a random optimization approach (Matya s, 1965; Sarma, 1990) (i) Choose a initial guess (b h , b w ) (ii) While some convergence criteria is not reached Of course, parameter determination could be pursued using more sophisticated methods, e.g. (Albani et al., 2021a) . In what follows we show the results coming out of parameter estimations. We then show results addressing how well the model predicts the number of new cases. We compare with real daily data for new cases and deaths, 3 and, to smooth out oscillations, we present the data using a 7-day moving average. One of the hurdles for compartmental models is the definition of initial conditions. Indeed, both the exact "first day" of the disease and how many people were infected in that first day are unknown. As a initial condition of our model, we assume that the first infection occurred 30 days before the official recordings. Indeed, the first day of the official data logs 489 cases; cf. (Lourenço et al., 2020) , who speculates that the epidemic in Italy and UK started one month before the first reported death. We also stipulate that there were 10 infectious at each age group at day one. Such uncertainty of initial conditions makes the results unreliable at the initial periods of simulations. In particular, the computation of (12) is not feasible for the lack of data. Thus, we impose m(j) ¼ 1 at the initial stages. We consider results for the City of Rio de Janeiro, henceforward simply called Rio. The results presented are for the sum of all ages: We recall the assumption that the contacts might lead to infections at home are time independent, and the contacts at work and other locations are identical, i.e., b w ¼ b o . Also, schools are closed, thus b s ¼ 0. We compare in Figs. 2 and 3 the results of our in-sample approximation of the dynamics of the disease. Note that the actual daily numbers of new cases (Fig. 2) and deaths (Fig. 3) are well approximated by the model. In Fig. 4 we display the values of b w and m, estimated by our random search. We see in particular that b w oscillates following the "macroscopic" oscillations of daily cases. Regarding the lethality strength m, note that it seems to decrease after day 90, in , except for the first one (that has a value roughly ten times bigger than the others). For visualization purposes we interpolated the discontinuous function b w using cubic splines. On the right we display the values of the lethality strength m, defined in (12). The blue curve is the least square quadratic curve. line with what was presented by (Dennis, McGovern, Vollmer, & Mateen, 2020; Horwitz et al., 2020) . Such aspect is captured by the concave least square parabola. Two important pieces of information coming out of the model are the values of the reproduction and replacement numbers R 0 and R t . Their computations follow easily from (4) and (5), and we display the results in Fig. 5 . Note that values of R 0 is characterized only by the population's profile and the contact matrix C e . On the other hand, R t depends also on the dynamics of the disease through the susceptible population. That is important since even if R 0 > 1, having R t smaller than one is enough to have a declining number of infectious. Another issue that is worth discussing is that we fit the number of clinically infective patients I c to the data. However, we discussed nothing thus far about the sub-clinical patients I sc . It turns out that the number of sub-clinical patients is, in this case, roughly ten times the number of registered patients, as shown in Fig. 6 . Such difference resonates well with claims that there is a significant sub-notification. For instance (Havers et al., 2020) , reports that the number of cases in the US might be at least 10 times of what is registered. It would not come as a surprise if this difference is even greater in Rio. Also, the WHO states that 80% of infections are mild or asymptomatic, 4 in particular among children, as discussed by (DeBiasi & Delaney, 2770; Han et al., 2019) . Recall that children are often not tested. Next, we investigate numerically if this model can make reasonable predictions. To do so, after the day 200, we run the simulations based on an in-sample (days 1e200) calculated b w and m, and compare with available data. We perform experiments considering the best and worst scenarios, i.e. lowest and highest number of cases and deaths. We also consider an intermediary case. The main idea shares similarities with some Value at Risk analyses used in risk management, where the past dynamics of investments help in predicting the future; see (Jorion, 2006) . We make the predictions by going back j 10-day periods and performing simulations with the highest and lowest values of b w and m in the period 200 À 10j and 200. Choosing j ¼ 4, we show in Figs. 7e9 the predicted number of daily and accumulated cases and deaths, with 70 days in advance. Note that, for most of the period of interest, the data is in between the best/worst case scenarios. We also plot the intermediary case, corresponding to the choice of b w and m being the average of the extreme values, often yielding more accurate predictions. We explore in this section simple immunization policies, shedding some light on the following question: if an age-based vaccination strategy should be implemented, who should be protected first? Certainly, real-life policies are never as simple as that, since other aspects as economics, logistics, ethics, politics, etc, play significant roles. 5 Not only that, but also it is not even clear that vaccinations schemes should be devised based on age groups (Pastor-Satorras & Vespignani, 2002) . However, vaccination policies based on age are easy to implement in practice. Based on the age-dependent reproduction and replacement numbers R 0;[ and R t;[ , we find out which age group contaminates most people. Consider for instance the dynamics of the disease in Rio, with a population profile as in Fig. 10 . Based on the evolution of the COVID-19 in the first 200 days in Rio, we computed the age-dependent replacement numbers R t;[ , displaying the results in Fig. 11 . Note that as far as spreading the disease is concerned, individuals of younger age are the most dangerous. And individuals above 60, which are the ones most affected by the disease, contribute less for its spreading. Fig. 6 . Number of sub-clinical cases in solid blue, corresponding roughly to ten times the number of official cases in dashed red. 4 https://www.who.int/emergencies/diseases/novel-coronavirus-2019/question-and-answers-hub/q-a-detail/q-a-similarities-and-differences-covid-19and-influenza. 5 For instance, free market policies are discussed in https://johnhcochrane.blogspot.com/2020/12/free-market-vaccines.html. Such heterogeneity leads to interesting immunization policy possibilities. People in the age group 15e19 (around 500 k people in Rio) have the most impact on virus dissemination, but few of them die. On the other hand, those above 60 (around 900 k people) are at high risk, but they have little impact on the spreading of the disease. To offer a glimpse on these issues, we modify our effective contact matrix C e , allowing for the possibility of immunization. Let us redefine (2) by where the transient diagonal matrix V(,) models the vaccinated age groups. Its j-th diagonal element is zero if the corresponding age group is protected, and one otherwise. Of course, it is also possible to consider partial immunization by setting values between zero and one. Consider the hypothetical situation of having b h and b w as computed before (see Fig. 4 ) and the population of Rio. Between the days 200e500, we assume that b w ¼ b opt ¼ 0.0064 or b w ¼ b max ¼ 0.0117, and a constant lethality strength m ¼ 1. Suppose that there is a limited supply of vaccines, which are made available only at a certain moment, that is, the vaccination campaign takes place at a single day. In our examples, this is 100, 150 or 200 days after the pandemic starts. Of course, the timing of vaccination campaigns impacts its effectiveness, and there is an important body of literature investigating this issue (Albani et al., 2021b; Amaku et al., 2021a; Quach and Deeks, 2021) . To better compare experiments, we define an utility function that takes into account the number of lost lives. For Y 2 (0, 1), we define the Y-discounted cost of death as X þ∞ t¼1 Y t dðtÞ; where d(,) is the number of daily fatalities. We consider Y ¼ 0.996 1/365 in the numerical tests. In our numerical experiments, we first assume that a single shot of the vaccine induces a perfect immunization and stop the virus spreading. In a second example, we consider partial immunization. Consider three age groups of similar size as targets of an immunization campaign: those aged 15e34 (z1,584 people), 40e59 (z1,636 people) and those above 50 (z1,670 people). In all tested situations, the groups with best results were either the 15e34 or those above 50. In Fig. 12 , we plot the accumulated deaths considering full immunization of the groups aged 15e34 or over 50, for a fixed b w ¼ 0.0064 after day 200. In solid blue we plot the total deaths if no vaccination takes place. The remaining plots correspond to the total deaths if vaccination occur on day 100 (black), 150 (red) or 200 (magenta). The solid lines correspond to vaccinating the young and the dashes lines corresponding to vaccination of the elderly. In Fig. 13 , we repeat the computations fixing a higher b w ¼ 0.0117 after day 200. Basically, two different policies are possible: stop the spreading corresponds to vaccinating those in the 15e34 age group, and protect the vulnerable corresponds to vaccinating those above 50. Note from the results that the impact of stop the spreading decreases as the vaccination date is delayed. That is because the percentage of young population already contaminated at later dates is larger, and vaccinating the young group ends up not being so effective. Such effect is amplified for higher transmission rates, as one can see by comparing Figs. 12 and 13. Note also that for short terms, it is always preferable to protect the vulnerable, indicating that vaccinating the elderly population might be a wise option if a continuous supply of vaccines is available. Considering the total amount of casualties in the long run, for lower transmission rates it is always better to vaccinate the elders (Fig. 12) . For higher transmission rates, these results are reversed (Fig. 13) . However, when considering the cost of death (14), the results are mixed, as shown in Table 1 . For related results, see (Meehan et al., 2020; 1093Giubilini, Savulescu, & Wilkinson, ; Zhao et al., 2020) . We also remark that although most countries are pursuing policies that prioritize the elderly, there are exceptions (The Japan Times, 2020). We consider in a second example the possibility that the vaccines are not perfect, and display the results in Fig. 14. Assume that a immunization campaign happens at day 150, and that either those aged 15e34 (solid lines) or those over 50 (dashed Fig. 8 . Similarly to Fig. 7 , the figure displays the number of daily deaths, dashed red being data and solid black being the best/intermediary/worst cases predicted by the model. Fig. 9 . We present here the accumulated number of cases (left) and deaths (right). In both figures, the dashed red plots correspond to real data and the solid black plots represent the simulation results for frozen b w and m, for the best/intermediary/worst case scenarios. plots) are vaccinated. Assume also that the immunization is perfect (in blue), or has partial efficiency of 70% (black) or 50% (red). Note that, as expected, as the vaccine efficiency decreases the number of casualties increase. It is maybe more surprising to realize that this effect is more aggravated for the elderly. The vaccination takes place on day 100 (black), 150 (red) or 200 (magenta). Note that vaccinating the elderly population is always optimal in the short term, but not necessarily in the long run. Fig. 13 . Accumulated deaths for b w ¼ 0.0117, a higher rate of transmission. As before, we assume no vaccination (solid blue), vaccination for group aged 15e34 (solid plots) and over 50 (dashed plots). The vaccination takes place on day 100 (black), 200 (red) or 300 (magenta). Note that vaccinating the elderly population is always optimal in the short term, but not necessarily in the long run. Epidemics as the COVID-19 are hard to model, not only because of biological factors, but also due to the unpredictability of human responses. Political, economical, social and individual factors influence the behavior of people, leading to various degrees of risk behavior. Making long term predictions is riskier than in general (non-human) biological systems, but they are crucial for planning non-pharmaceutical interventions, allocation of resources, economical decisions, etc. We propose in the work a general form to predict possible values of crucial parameters based on their history. Our method has two steps. Based on the dynamics of the disease on a multi-generational SIR model, we first determine the values of some time-dependent parameters, using a random optimization algorithm. Then, we feed the SIR model using extreme values of such parameters, and predict best/worst case scenarios. We are also consider a prediction that is in between the two extreme cases. Our methodology is general and can be applied to other systems and circumstances, and might be useful to make longterm predictions in situations where unexpected occurrences might change the behavior of the system. We also explore age-based immunization strategies, under simplifying hypotheses. For instance, while we assume the behavior of the population as known, the very existence of vaccines changes people's behavior. Moreover, several important aspects of vaccination are not considered here, as ethical problems, logistics, etc. We also point out that other policies might be better, as targeting certain classes of workers or individuals with comorbidities. In our limited setting, we find out that optimal policies depend, among other things, on transmission rates, the vaccination efficiency and its timing. Earlier vaccination dates, low vaccine efficiency and higher transmission indices increase the impact of vaccinating the young. Under most conditions however, it pays off to vaccinate the elderly, although they contribute the least for spreading the viruses. Indeed, for short terms, protecting the older population is always preferred. That seems to indicate that under a continuous supply of vaccines, protecting the elderly is a wise option. It is clear that the model provides different indications under different scenarios, and we shall explore more detailed policies in a future work. As we point out in our Disclaimer, our goal is not to predict actual number of cases and deaths or to propose general policies, but instead to investigate limited aspects of the disease dynamics and raise possibilities of immunization strategies. Table 1 Numerical investigation of immunization policies based on age. In each line we present the values of the death cost functional (14) assuming different age groups are vaccinated. The second column contains the population of the groups. In the three columns that follow, we present the results considering the vaccination takes place on the days 100, 150 or 200 and assuming that b w ¼ 0.0064. The following three columns display the results for b w ¼ 0.0117 and vaccination on the days 100, 150 or 200. For each column, we highlight the best results. The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. The interplay between voluntary vaccination and reduction of risky behavior: A general behavior-implicit sir model for vaccine preventable infections A multi-risk sir model with optimally targeted lockdown Herd immunity thresholds for sars-cov-2 estimated from unfolding epidemics The impact of covid-19 vaccination delay: A modelling study for Chicago and NYC data Estimating, monitoring, and forecasting covid-19 epidemics: A spatiotemporal approach applied to nyc data A simple planning problem for COVID-19 lockdown Modelling the impact of delaying vaccination against sars-cov-2 assuming unlimited vaccines supply Modelling the test, trace and quarantine strategy to control the covid-19 epidemic in the state of São Paulo, Brazil. Infectious Disease Modelling, 6, 46e55 Towards providing effective data-driven responses to predict the Covid-19 in São Paulo and Brazil What will be the economic impact of COVID-19 in the US? Rough estimates of disease scenarios, working paper 26867 Epidemiological short-term forecasting with model reduction of parametric compartmental models. application to the first pandemic wave of COVID-19 in France On the basic reproduction number in continuously structured populations COVID-19 modelling: The effects of social distancing Hyperbolic models for the spread of epidemics on networks: Kinetic description and numerical methods Sex differential in covid-19 mortality varies markedly by age A macroeconomia das epidemias: Heterogeneidade interestadual no brasil Mathematical epidemiology: Past, present, and future COVID-19 infection fatality ratio: Estimates from seroprevalence A mathematical model reveals the influence of population heterogeneity on herd immunity to sars-cov-2 An economic model of the COVID-19 epidemic: The importance of testing and age-specific policies Estimating the extent of asymptomatic covid-19 and its potential for community transmission: Systematic review and meta-analysis Bayesian dynamical estimation of the parameters of an SE(A)IR COVID-19 spread model The economics of vaccination Modeling the post-containment elimination of transmission of COVID-19 Influence of non-homogeneous mixing on final epidemic size in a meta-population model Symptomatic and asymptomatic viral shedding in pediatric patients infected with severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2): Under the surface Complexity of the basic reproduction number (R 0 ) Improving COVID-19 critical care mortality over time in England: A national cohort study On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations Modeling the heterogeneity in COVID-19's reproductive number and its impact on predictive scenarios The macroeconomics of epidemics, Working Paper 26882 Impact of non-pharmaceutical interventions (NPIs) to reduce COVID-19 mortality and healthcare demand Nine challenges in incorporating the dynamics of behaviour in infectious diseases models COVID-19 vaccine: Vaccinate the young to protect the old? Individual variation in susceptibility or exposure to sars-cov-2 lowers the herd immunity threshold On the spatial dynamics of vaccination: A spatial SIRS-V model Macroeconomic implications of COVID-19: Can negative supply shocks cause demand shortages? The milroy lectures on epidemic disease in englandethe evidence of variability and of persistency of type Clinical characteristics and viral RNA detection in children with coronavirus disease 2019 in the Republic of Korea Seroprevalence of antibodies to SARS-CoV-2 in 10 sites in the United States COVID-19: What proportion are asymptomatic The mathematics of infectious diseases Trends in covid-19 risk-adjusted mortality rates in a single health system The basic approach to age-structured population dynamics Modeling the control of COVID-19: Impact of policy interventions and meteorological factors Gender differences in patients with COVID-19: Focus on severity and mortality. Frontiers in Public Health Value at risk: The new benchmark for managing financial risk A contribution to the mathematical theory of epidemics Age-structured model for COVID-19: Effectiveness of social distancing and contact reduction in Kenya The stage-structured epidemic: Linking disease and demography with a multi-state matrix approach model Stage-structured transmission of phocine distemper virus in the Dutch Simulando os efeitos de políticas de distanciamento social em São Paulo usando um modelo SEIR Integrating behavioral choice into epidemiological models of AIDS Global dynamics of an SIR epidemic model with nonlocal diffusion Analytic models for SIR disease spread on random spatial networks Assessing the age specificity of infection fatality rates for covid-19: Systematic review, meta-analysis, and public policy implications Linear age-structured population models as a base of age-structured epidemic models Fundamental principles of epidemic spread highlight the immediate need for large-scale serological surveys to assess the stage of the sars-cov-2 epidemic Modeling the interplay between human behavior and the spread of infectious diseases Random optimization Spatial heterogeneity and the design of immunization programs Age-targeted dose allocation can halve covid-19 vaccine requirements Applications of mathematics to medical problems Social contacts and mixing patterns relevant to the spread of infectious diseases Continuous and discrete sir-models with spatial distributions Immunization of complex networks Modeling future spread of infections via mobile geolocation data and population dynamics. an application to COVID-19 in Brazil Optimal control of the COVID-19 pandemic with non-pharmaceutical interventions Merging economics and epidemiology to improve the prediction and management of infectious disease Projecting social contact matrices in 152 countries using contact surveys and demographic data The effect of control strategies to reduce social mixing on outcomes of the COVID-19 epidemic in Wuhan, China: A modelling study Covid-19 vaccination: why extend the interval between doses? Why is it difficult to accurately predict the COVID-19 epidemic? Infectious Disease Modelling, 5, 271e281 On the convergence of the Baba and Dorea random optimization methods Using insights from behavioral economics to mitigate the spread of COVID-19. Applied Health Economics and Health Policy Global stability of the endemic equilibrium of multigroup sir models with nonlinear incidence High order discretizations for spatial dependent SIR models With goal of herd immunity, Indonesia will give younger people vaccines first Social contact patterns and control strategies for influenza in the elderly Estimates of the severity of coronavirus disease 2019: A model-based analysis. The Lancet Infectious Diseases Simulating the spread of COVID-19 via spatially-resolved susceptibleexposed-infected-recovered-deceased (SEIRD) model with heterogeneous diffusion The impact of COVID-19 and strategies for mitigation and suppression in low-and middle-income countries Sars-cov-2 detection, viral load and infectivity over the course of an infection Modeling and forecasting COVID-19. Notices of the Asymptotic profiles of the endemic equilibrium to a diffusive sis epidemic model with mass action infection mechanism Characteristics of and important lessons from the coronavirus disease 2019 (COVID-19) outbreak in China: Summary of a report of 72314 cases from the Chinese center for disease control and prevention Asymptotic profiles of steady states for a diffusive SIS epidemic model with mass action infection mechanism The optimal vaccination strategy to control COVID-19: A modeling study based on the transmission scenario in Wuhan city Duration of viral shedding of discharged patients with severe COVID-19 The authors would like to acknowledge the financial support of CNPq and FAPERJ. This study was financed in part by the Coordenação de Aperfeiçoamento de Pessoal de Nível Superior -Brasil (CAPES) -Finance Code 001.