The next Zika virus outbreak will affect young women of reproductive age disproportionately: mathematical modeling study

9 Zika virus (ZIKV) is a mosquito-borne flavivirus. In 2016, Managua, the capital city of 10 Nicaragua experienced a Zika virus outbreak. Amongst confirmed cases, seropositivity increased 11 with age and >50% of the adult population had ZIKV antibodies. The duration of immunity 12 to ZIKV remains unclear. We used an agent-based model in a Bayesian framework to describe 13 the dynamics of immunity to ZIKV, and predict the future age-specific risk of ZIKV infection 14 in the population of Managua. We also investigated the potential impact of a ZIKV vaccine 15 under different assumptions about the persistence of immunity. Assuming lifelong immunity, 16 the risk of a ZIKV outbreak will remain low until 2035 and rise above 50% in 2047. Young 17 women of reproductive age will be at highest risk of infection during the next ZIKV outbreak. 18 ZIKV vaccine development and licensure are urgent, to attain the maximum benefit in reducing 19 the population-level risk of infection and the risk of adverse congenital outcomes. 20

inferences from data about the 2016 ZIKV epidemic (Table 1). We considered different scenarios 87 about the duration of immunity, the timing and scale of ZIKV reintroductions in the population, 88 and the timing and scale of a hypothetical vaccination program targeted towards 15 year old girls. 89 2.2 Model structure 90 We simulated a population of 10,000 individuals for 80 years . We assigned agents age 91 and ZIKV infection status (susceptible S, infected I or immune R). Initial conditions reflect the 92 situation in Managua, Nicaragua in 2017, when there was no documentation of active transmission.

93
In the outbreak-free period, we only considered demographic and immunity processes: births, deaths, 94 ageing and, if applicable, loss of immunity and vaccination. Given the scarcity of these events 95 at the individual level, we select a long time-step of seven days and stochastically evaluated the 96 epidemic-related transition probabilities at each time step for each agent. After a given time, ZIKV 97 cases were reintroduced in the population. Upon reintroduction, the time step was reduced to 0.1 98 days, and we evaluated the transition probabilities. Susceptible agents may become infected at 99 a rate β a I/N , where β a is the age-dependent transmission rate and N the total population size. 100 We ignored the vector population and assumed that the force of infection is proportional to the 101 overall infection prevalence. We allowed six months for the outbreak to finish after introduction.
times. In the baseline scenario, we assumed no vaccination, no loss of immunity and an introduction of 10 infected individuals per simulation.
We inferred the probability distributions for the age-specific transmission rate β a and the recovery

126
We ignored demography in this model because it covers a short time span. We recorded the overall 127 cumulative incidence of ZIKV cases using a dummy compartment: in order to compute the weekly incidence on week t: We fitted the model to weekly incidence data A using a normal likelihood after a square-root 132 variance-stabilizing transformation (Guan, 2009): where ρ is a reporting rate parameter and σ an error parameter. In addition, we also fitted the 135 model to the number of individuals with anti-ZIKV antibodies at the end of the epidemic by age 136 group B a using a binomial likelihood: 138 where B a the number of individuals with antibodies, n a is the number of individuals in each age 139 group, and p a = R a (t end )/N a (t end ) the proportion of immune at the end of the epidemic. The 140 full likelihood was obtained by multiplying Eq. 6 and Eq. 7. We chose weakly-informative priors for all parameters and fitted the model in Stan (Table 2). We describe the calculation of the 142 basic reproduction number R 0 in appendix A.2. We used one thousand posterior samples for β a   148 We used the deterministic model, described in the previous section, to infer the proportion of people

153
Given the absence of evidence about the duration of immunity to ZIKV, we considered four scenarios 154 assuming exponentially distributed durations of immunity with means of 30, 60, 90, or 150 years.

155
These values correspond to a proportion of initially immune agents that loses immunity after 10 156 years of 28%, 15%, 11% or 6%, respectively (Appendix A.3).   165 We reintroduced ZIKV in the population after a delay of d = {1, · · · , 80} years in independent 166 simulations. We chose this approach rather than continuous reintroductions to remove some of 167 the stochasticity and assess more clearly the association between immunity decay and risk of an

193
From the simulations, we collected 1) the evolution of the age-specific ZIKV immunity in the pop-

208
In our simulation, the expected population size increased by 42% between 2017 and 2097. Under 209 the assumption that ZIKV infection results in lifelong protective immunity, population renewal 210 will create an imbalance in the proportion immune in different age groups. We expect the overall 211 proportion of protective immunity in the population to have halved (from 48% to 24%) by 2051 and 212 to be concentrated among the older age classes ( Fig. 2A). The 0-14 year old age group will become 213 entirely susceptible by 2031 and the 15-29 year old age group by 2046. 2B). After this point, attack rates resulting from ZIKV reintroduction will rise steeply. By 2047, we 218 predict that ZIKV reintroductions will have a 50% probability of resulting in outbreaks with attack 219 rates greater than 1% (Fig. 2C).

Loss of immunity 234
If protective immunity to ZIKV is not lifelong, the window before observing a rise in the attack rates 235 resulting from ZIKV reintroduction will shorten (Fig. 4A). For instance, if 15% of the those who 236 were infected in 2016 lose their immunity after 10 years (a mean duration of immunity of 60 years), Nicaragua, after reintroduction, will remain low before 2035 because of herd immunity, then rise to 254 50% in 2047. If protective immunity to ZIKV decays with time, which remains to be elucidated,

255
ZIKV recurrence could occur sooner; in that case, the distribution of the risk would be more equal   age distribution of future ZIKV outbreaks will likely differ and that a new large epidemic will be 309 delayed for "at least a decade", (Ferguson et al., 2016). Our work, which models the processes of 310 ageing in relation to protective immunity to ZIKV explicitly, provides a quantified confirmation of 311 these previous works, but adds the effect of different transmission rates by age group.

342
The generalizability of some conclusions is limited. ZIKV has touched most tropical regions of 343 the Americas. In many of these regions, ZIKV epidemics resulted in a high post-epidemic seropreva-344 lence, similar to Nicaragua. However, differences in vector density and demography limit the scope 345 and generalizability of our predictions regarding the time until ZIKV recurrence. Conversely, the 346 increased burden experienced by younger age groups and the subsequent higher-than-expected num-347 ber of adverse congenital events resulting from a ZIKV outbreak is likely to hold in many settings.

348
Seroprevalence studies across the Americas in all age strata, using tests with high specificity, would 349 greatly help risk prediction.

351
Preparedness is vital; the time until the next outbreak gives us to opportunity to be prepared.

352
The next sizeable ZIKV outbreak in Nicaragua will likely not occur before 2040 but the probability 353 of outbreaks will increase. Young women of reproductive age will be at highest risk of infection Here, we provide the pseudo code of the ABM (Algorithm 1).

A.2 R 0 522
We used the next generation matrix method described by Diekmann et al. to calculate R 0 (eq. 8 523 -10). β 1 is the transmission rate for the 0-14 age group; β 2 for the >15 group; γ is the common 524 recovery rate.
A.3 Loss of immunity scenarios 526 We explored plausible scenarios of loss of immunity ( Fig. A.1). A.6 Comparison of SIR model with SEIR model and the Pandey model. 541 We compared the SIR model with a SEIR and model that explicitely models the vector; the Pandey of the more complex models does not outperform the fit of the simplest (SIR) model (