key: cord-0902044-p70a6qg1 authors: Chen, Chun-Miin (Jimmy); Stanciu, Alia C. title: Simulating influenza epidemics with waning vaccine immunity date: 2021-09-10 journal: Medicine (Baltimore) DOI: 10.1097/md.0000000000027169 sha: 23778261a397eb15a7996fefffb798dd415c8d1b doc_id: 902044 cord_uid: p70a6qg1 Observational studies indicate that vaccine-induced immunity can decline over time. However, few researchers have incorporated this kind of waning effect into their virus spread models. In this study, we simulate an influenza epidemic that considers the effects of waning immunity by fitting epidemiological models to CDC secondary historical data aggregated on a weekly basis, and derive the transmission rates at which susceptible individuals become infected over the course of the influenza season. Using a system of differential equations, we define four groups of individuals in a population: susceptible, vaccinated, infected, and recovered. We show that a larger number of initially infected individuals might not only bring the influenza season to an end sooner but also reduce the epidemic size. Moreover, any influenza virus that entails a faster recovery rate does not necessarily lead to a smaller epidemic size. We illustrate how simulation helps in understanding the effects of influenza epidemiological model in the presence of waning influenza vaccine immunity. Influenza is a contagious disease around the world. Year after year, influenza viruses cause illness to people of all ages, which can lead to influenza-related complications or even trigger more severe inflammation in the body. In fact, influenza has been ranked as the eighth leading cause of death in 2017. [1] According to the World Health Organization, [2] annual influenza-related death tolls ranged from 290,000 to 650,000 worldwide. Since 1997, the Centers for Disease Control and Prevention (CDC) have been tracking influenza activity and influenza-related illnesses. In the United States, CDC produces influenza surveillance data, such as the percentage of specimens testing positive for an influenza virus. Every week, all clinical laboratories report to CDC the total number of respiratory specimens tested and the number that tested positive for influenza viruses. Recently, researchers have started to thoroughly investigate and give special attention to the waning effectiveness of vaccination-induced immunity in the human body. In the 2017 to 2018 influenza season, there were as many as 50 viruses identified for vaccine productions. [3] Given the selected viruses for the new formulation, vaccine manufacturers are racing against time as they produce, test, and distribute the vaccine under a tight schedule. [3] [4] [5] Unfortunately, antigenically variable pathogens, such as influenza viruses, are living organisms and are capable of continuously altering the proteins or carbohydrates on their surface and mutating to different strains, thereby potentially escaping from the immunity induced by vaccination. [6] Many studies have reported that vaccine effectiveness within a single influenza season can decrease over time. Kissling et al [7] stated that the decreasing influenza vaccine effectiveness could be the result of virus changes through the season or waning immunity. Pebody et al [8] and Belongia et al [9] advocated the need for developing influenza vaccines that provide better and longerlasting protection in view of the occurrence of late past-season outbreaks. Hill et al [10] developed a multi-strain transmission model to predict the chronological interactions between the influenza viruses and waning immunity among the population. In this study, we aim to update classic epidemiological models by explicitly accounting for waning vaccine immunity over time, as well as to examine the effects that various parameters have on the size of the influenza epidemic. Methodologically, we employ simulation, which has been one common approach for understanding the course of an influenza epidemic both in practice and in the literature. [11, 12] Goeyvaerts et al [13] developed a seasonal influenza transmission model for evaluating the impact of vaccination on the incidence of infection, disease, and mortality. The authors proposed to directly estimate the values of the seasonal influenza model parameters by fitting the model to influenza data over multiple influenza seasons. To empirically ground the model, we follow the literature and present a simplified, yet accurate way of fitting influenza seasonal transmission models to publicly available data. [13, 14] Using a system of differential equations, we simulate the dynamics of the infection in the population compartmentalized into several groups. The differential equations allow us to define the rates at which susceptible individuals become infectious or infected individuals recover from the disease. Lastly, we conduct an extensive sensitivity analysis to validate our model by varying a wide range of common experimental parameters discussed in the extant literature. This study visualizes the progress of the influenza epidemic not only when vaccine effectiveness wanes, but also when some seemingly favorable change of conditions counter-intuitively contributes to a worsening of the spread of influenza. The remainder of this article is organized as follows: Section 2 presents the classic susceptible-infectious-recovered (SIR) model that is often used for segmenting a population into groups during an influenza epidemic. In Section 3, we conduct simulation studies using various scenarios for investigating the implications of the empirically based periodic transmission rate, waning vaccination effectiveness, and viciousness of the influenza virus on the size of the epidemic. Finally, in Section 4, we conclude with a summary of our findings. In this section, we present an approach for fitting a model to influenza data. First, we define the groups that each individual in a population will belong to during an influenza epidemic. Then, using differential equations, we devise a deterministic epidemic model to investigate the implications of varying model parameter values on the severity of the influenza epidemic in the population. Using a sinusoidal function for simulating influenza transmission is not uncommon in epidemiological research. [15] [16] [17] For a given empirical data set, we fit the following sinusoidal function to the weekly percentage of the positive laboratory specimens test results over the course of a year: where A specifies the amplitude of the oscillation, B defines the time horizon of a complete season, C determines the starting week in a season, D shifts f(t) such that it is positive for all integers t∈ 1; 2; . . . ; 51; 52 f g , and t indexes the week number in a given year. We then construct the following optimization problem to get the fitted parameter values for modeling the seasonal influenza phenomenon: We use x(t) to denote the weekly percentage of specimens with positive results in any available historical data set. Once we obtain the optimal parameter values for the decision variables A, B, C, and D, we can use the optimized sinusoidal to derive a contact function that will simulate the rate at which an infected individual transmits the disease to other susceptible individuals. The contact function shows how soon individuals in the population can move between the susceptible, infectious, and recovered groups. In particular, the contact (transmission) rate at which a susceptible individual becomes infected can be obtained by taking the first derivative of Eq. (1) with respect to t. The contact rate, denoted by bðtÞ ¼ f 0 ðtÞ, is a cosinusoidal function of time that informs the various rates of transmission during the influenza season. Many studies have used differential equation-based simulation models for simulating the course of influenza. [18] [19] [20] Typically, these models place individuals in a predetermined number of groups. In this study, we assume that, at any given time, each individual in the population (N) belongs to one (and only one) of the following groups: susceptible (S), vaccinated (V), infectious (I), or recovered (R). [11, [21] [22] [23] [24] In particular, S includes those who are healthy but exposed to the risk of contracting influenza from any infectious individuals; V refers to the proportion of individuals in N who are being vaccinated but who can also continuously leak to group S; I denotes those who are infectious agents of influenza to those in S; and R denotes those who have contracted influenza, recovered, and are therefore immune for the remainder of the season. Figure 1 exhibits a flowchart that conceptually illustrates to which group the individuals in a population could belong during an influenza epidemic. Note that the dotted line between V and S is the connection that has not yet been considered in most of the influenza-related research so far. Given the definition of the SVIR epidemic model, we developed a system of deterministic differential equations that characterizes the rates at which any individual in one group can move to another group [25] : We use z to denotes the rate of the decreasing vaccine effectiveness (i.e., waning effect). Studies showed that the vaccineinduced antibody titers can decrease, within a year after the inoculation, to levels encountered in those without protection. [26] Chen and Stanciu Medicine (2021) 100: 36 Medicine Another biological definition of the waning effect is that the antigens on the influenza virion surface (i.e., the target of protective immunity) continually evolves and mutates in the human population. Hence, the influenza strain increasingly differs from its progenitor so that the vaccine-induced immunity against the progenitor can diminish over time. [27] Mathematically, we assume some percent of the vaccinated group can lose the immunity protection and rejoin to the susceptible group. Additionally, g denotes the inverse of the average length of time an individual stays infectious until moving to the R group. For example, if infected individuals take, on average, four days to recover, then g ¼ 0:25. To solve the differential equations, we need the following initial conditions: where t 0 denotes the time when the pathogen is introduced in the population, and a represents the level of population immunity given the vaccine's effectiveness and the proportion of susceptible individuals being vaccinated. Admittedly, some more advanced, complex influenza transmission models with additional parameters will be able to better capture other temporal or systematic variations in the time series data. For example, while the differential equations are popular in simulating the course of influenza transmission, they tend to focus more on modeling the behavior of the groups as a whole. To focus on modeling the probability of individuals being infected given their demographic or community information, agent-based models are capable of considering the stochastic nature of the transmission rate between individuals. [11, 24] For more discussions on the differential equations and the agent-based simulation models, see Paleshi et al [12] Our parsimonious mathematical models and differential equations are meant for capturing the main patterns that are evident in the empirical data for the following numerical study. From the various publicly available data sets on the CDC's website, we focus on the weekly percentage of the positive laboratory specimens test results and use those as a proxy for the fraction of the population that is infected. The data were collected from ∼300 clinical laboratories located throughout all 50 states, Puerto Rico, Guam, and the District of Columbia. Prior to the 2015 to 2016 influenza season, the weekly influenza update from public health and clinical laboratories were not separately reported. Compared to the clinical laboratories, public health laboratories often receive samples that have already tested positive for an influenza virus at a clinical laboratory. [3] Relatively speaking, the data from clinical laboratories are less distorted than those from public health laboratories in terms of truthfully representing the situations in the population. Therefore, we chose to use data from four influenza seasons, starting with the 2015 to 2016 season, to inform the estimates for our epidemiological model. Figure 2 displays the time-series data of four (over-imposed) influenza seasons consisting of 208 (52 Â 4) weekly fractions of the positive laboratory specimens test results. Overall, the percentages rise rather quickly at the beginning of the season before declining to some low levels starting with week 24. Note that every plot has peaks between week 52 and week 12 and a small uptick before the end of the season, at around week 36. Because we use only deidentified statistical data aggregated by the CDC at the national level on a weekly basis, no information can be used to uniquely identify each individual. Thus, this study does not involve human subjects in a way that would require review or approval by Bucknell University's Institutional Review Board. As opposed to directly imposing some function for the influenza contact (transmission) rate, [20] we derive the cosinusoidal contact rate from the empirical data. To this end, we solved the minimization problem in Eq. (2) using MATLAB to obtain A, B, C, and D and specify Eq. (1) (solver sequential quadratic programming; www.mathworks.com). In the definition of x(t) in Eq. (2), we use the average of the percentages in Figure 2 , that is, 52 weekly average ratios across the four influenza seasons. (e.g., x (1) is the average of the four first-week percentages from the data sets.) For the given empirical data set, we fit the sinusoidal function (1), using the model specified in (2) Next, we derive the estimated contact rate function at which a susceptible individual becomes infected when contacting an infectious individual. Intuitively, the contact rate determines the instantaneous rate of change to the weekly positivity rate. We propose a novel approach for finding the contact rate by borrowing a fundamental concept from Physics: the contact rate is to the weekly positivity rate as the acceleration is to velocity. That is, if we view the positivity rate in a week as velocity, then the contact rate should be the instantaneous rate of acceleration in the given week. Thus, we take the first derivative of the sine function f(t) with respect to t to find the contact rate as a cosine function: bðtÞ ¼ 1:3764 cos 2p 43:0102 ðt À 139:2729Þ þ 1 As such, the contact function Eq. (7) has a maximum at week 2 (in January) and a minimum at week 24 (in June) for a given year, matching the empirical observation that influenza activity often begins in October, peaks sometime between December and February, and lasts as late as May. [3] Note that the constant 1 in Eq. (7) is added to keep the functionbðtÞ realistic for all integer t∈ 1; 2; . . . ; 51; 52 f g . Without the added constant, Eq. (7) oscillates between +1.3764 and À1.3764 over t. If adding any constant <1, then Eq. (7) might render some negative contact rate. If the added constant is >1, then Eq. (7) might render some increasing infection rate for the population even during the summertime. This study also examines the impact of the length of the infectious period on the size of the epidemic. According to CDC, [3] infected people can be contagious one day before experiencing symptoms, which typically last 3 to 7 days, and can stay contagious for up to 7 days after becoming sick. The contagious period can be longer than 7 days for children and some people with weakened immune systems. Cori et al [28] pointed out that it is impossible to directly observe the duration of the infectious period. The authors reviewed the literature and found that the 95% probability intervals of the mean duration of the infectious period can be any number of days up to 12 days, while the majority of cases having an infectious period are shorter than 2.9 days. Thus, we consider g, the recovery rate of infected individuals, to vary between 1 12 and 1 2 . Since it is of interest, we also tracked and observed during our simulations the total number of recovered individuals, as a proxy for the size of the epidemic. Another uncertainty that we bring into the study is the proportion of people who gain immunity, at least in the beginning of the season. Statistics from the 2018 influenza season showed that, depending on age group, the percentage of people who received an influenza vaccination varies between 30% and 70%, with vaccine effectiveness ranging from 10% to 60%. [29, 30] Thus, we let a, the proportion of susceptible individuals fully protected by the vaccine, vary between 3% and 42%. Moreover, we vary I 0 , the initial number of infected individuals, between 1 and 100. Lastly, Grohskopf et al [31] pointed out that the observed waning effect might not consistently manifest itself across all age groups, virus subtypes, or seasons. Nonetheless, Ferdinands et al [32] examined the association between influenza vaccine effectiveness and time since vaccination among patients in the United States. The authors found significant evidence that influenza vaccine protection could decline by 6% to 11% per month since the time of vaccination, depending on influenza type. Thus, we vary z, the rate at which the number of (newly) vaccinated individuals decreases over time, from 6% to 11%. We code the system of differential Eqs. (3) to (6) in MATLAB and solve them numerically for the values of a given t under different scenarios (solver ode45; www.mathworks.com). When computing, we begin by assigning the initial values of the groups of the population, with N and R 0 being always held constant at 10 6 and at 0, respectively. In addition, every influenza season spans 52 weeks, starting at week 40 of 1 year and ending at week 39 of the following year. During the season, individuals continuously interact with each other, and the size of the four groups changes at the rates defined by the differential equations. Figure 4 shows the impact of the recovery rate on the number of individuals in the four groups over the course of the influenza season. We observe that an influenza virus that requires a longer recovery time tends to cause the population to experience another uptick in the number of infected individuals at around week 36. We also note that the longer the time the infected individuals take to recover, the larger the epidemic size is-that is, the total number of recovered individuals at the end of the season. Figure 5 shows that if the population has a larger number of initially infected individuals (I), then the peak of the season would not only arrive sooner but also be smaller. A subtle and counterintuitive observation is that the epidemic size-that is, the total number of recovered individuals at the end of the season-is actually smaller when I is larger. We will examine this phenomenon in more detail shortly. Figure 6 shows the impact of the proportion of the vaccinated population on the number of individuals in the four groups over the course of the influenza season. The greater the proportion of individuals receiving the vaccination at the beginning of the season, the smaller the epidemic size (the total number of 40 52 12 24 36 Week of the year Chen and Stanciu Medicine (2021) 100: 36 Medicine recovered individuals) at the end of the season. An increase in a seems to lead to a decrease in the size of the epidemic, but it also delays the peak of the season, as well as leading to an uptick in infected individuals toward the end of the season. Figure 7 highlights the missing link in the literature-the effect of the individuals leaking from V to S when the vaccination is subject to the waning effect. Even with a slight increase in the waning effect, for example, z from 0 to 6% (i.e., top-left panel to top-right panel), the size of the epidemic could increase by more than 100%. Therefore, researchers studying related topics should not disregard the waning effect. Although a higher waning rate leads to a larger epidemic size, the peak of the season takes place roughly in the same week. The implication is that researchers or practitioners may not be able to tell if the effectiveness of a given vaccine can wane over time by observing the timing of the peak of the season. Figure 8 shows the epidemic size given four different values for the number of initially infected individuals. Counter-intuitively, relatively faster recovery rates, roughly between 1/12 and 1/7, can actually make the epidemic worse. Another interesting observa-tion is that, for the cases where the recovery rate is greater than about 15%, the greater the number of initially infected individuals, the smaller the epidemic (i.e., the smaller the total number of recovered individuals at the end of the season). Note that for influenza that entails a recovery period of 7 days or fewer, a greater number of initially infected individuals can lead to a decrease in the season's epidemic size. But for influenza that entails a recovery period of 7 days or more, a greater number of initially infected individuals can lead to an increase in the season's epidemic size. Figure 9 illustrates that the higher the waning rate, the larger the epidemic size. The implication is that the effectiveness of the vaccine must be made as good as possible to effectively mitigate the potential outbreak of an epidemic. Similar to the phenomenon observed in Figure 8 , a faster recovery rate might increase the epidemic size, and the phenomenon holds even under different rates of the waning effect. To conduct a sensitivity analysis, we separately fitted the models using the data of the individual seasons and obtained four different sets of values of the parameters A, B, C, and D as shown Table 1 . Not reported here, all of the results using the parameter values led to the same conclusions as do those using the average of the four seasons' data. Influenza has been a leading cause of death in the United States. Despite vaccination, the number of susceptible, infected, and recovered individuals seems to follow similar patterns year after year. [24] In this study, we combined multiple approaches to better capture the complexities of the influenza epidemic phenomenon. The theoretical contribution of our study lies in fitting empirical influenza data to the models that capture the periodic phenomena of influenza infections and the oscillating contact rate over time. Moreover, we extended the classic SIR model to the SVIR model by incorporating the effect of waning vaccination immunity to better reflect the course and spread of influenza infections. The 40 52 12 24 36 Week of the year Chen and Stanciu Medicine (2021) 100: 36 Medicine practical implication of the study lies in simulating the size of the influenza epidemic when varying several model parameters. In particular, we show how an uptick in influenza cases toward the end of the season or even a delay in the peak of the season might be caused by scenarios involving an influenza virus that entails a longer recovery time, or when a greater portion of the population gets vaccinated. The rather counter-intuitive results suggest that a larger number of initially infected individuals might not only bring the influenza season to an end sooner but also reduce the size of the epidemic. Moreover, an influenza virus that entails a faster recovery rate does not necessarily lead to a smaller epidemic size. One limitation of the study may be attributed to the availability of the less distorted data we chose to use, which was unavailable on the CDC website until 2015. In addition, the differential equation model assumes that all subjects are directly connected with each other, and as likely to get or spread the disease, which may not be true for populations that are sparsely connected. Nevertheless, the numerical experiments present a number of interesting phenomena that are worth pursuing further. Future research and extensions of the investigations in this study can consider different types of influenza viruses. According to CDC, [3] both type A and type B influenza viruses are responsible for causing seasonal epidemics of the disease every winter. [10] It has been reported that different types of the influenza virus could show up at different times in an influenza season. [33] Thus, the implications of the waning effect over time we discussed in this study may pan out differently for different types of the influenza virus. Another direction for future research is to consider the group of immunosuppressed individuals in the model. [34] Individuals with immune dysfunction may not respond to the vaccine as well as people without immune dysfunction. [35] Furthermore, people who are severely immunocompromised (e.g., HIV-infected people, organ transplant recipients, or people who undergo immunosuppressive treatments) tend to prolong influenza virus shedding, increase morbidity or mortality following infection, and be reinfected with the same virus strain. [36, 37] Thus, researchers can study the impact of immunosuppression on the SVIR model. Finally, researchers may consider applying the framework of this study and analyze influenza data from other countries or regions, for additional insights. As the entire world is increasingly inoculated with COVID-19 vaccines, we anticipate that many more studies will be underway to further model the immunity provided by the vaccines. The waning effect we consider in this article, even though introduced in the context of an influenza epidemic, could become especially relevant in modeling and understanding the COVID-19 pandemic immunity-related parameters. At the same time, the reliability and the frequency of data reporting will have a direct effect on gaining a better understanding of vaccine administration. Writingoriginal draft: Chun-Miin (Jimmy) Chen, Alia C. Stanciu. Writingreview & editing: Chun-Miin (Jimmy) Chen, Alia C. Stanciu. Deaths: Final Data for 2017 Ask the expert: Influenza Q&A Influenza (Flu) Global challenges in seasonal influenza vaccine supply, use, and policy Influenza vaccine response during the start of a pandemic: report of the Third WHO informal consultation Mapping the antigenic and genetic evolution of influenza virus Low and decreasing vaccine effectiveness against influenza A(H3) in 2011/12 among vaccination target groups in Europe: results from the I-MOVE multicentre casecontrol study Vaccine effectiveness of 2011/12 trivalent seasonal influenza vaccine in preventing laboratoryconfirmed influenza in primary care in the United Kingdom: evidence of waning intra-seasonal protection Waning vaccine protection against influenza A (H3N2) illness in children and older adults during a single season Seasonal influenza: modelling approaches to capture immunity propagation Simulation of strategies for containing pandemic influenza Simulation of mitigation strategies for a pandemic influenza Estimating dynamic transmission model parameters for seasonal influenza by fitting to age and season-specific influenza-like illness incidence Is the impact of childhood influenza vaccination less than expected: a transmission modelling study Modelling antigenic drift in weekly flu incidence Mortality due to influenza in the United States-an annualized regression approach using multiple-cause mortality data Estimating influenza deaths in Canada Antiviral treatment for pandemic influenza: assessing potential repercussions using a seasonally forced SIR model Pandemic H1N1 influenza: predicting the course of a pandemic and assessing the efficacy of the planned vaccination programme in the United States Dynamical resonance can account for seasonality of influenza epidemics A contribution to the mathematical theory of epidemics Accounting for individual behaviors in a pandemic disease spread model Implications of vaccination and waning immunity A primer for hybrid modeling and simulation Modeling the effects of vaccination and treatment on pandemic influenza Effectiveness, immunogenicity, and safety of influenza vaccines with MF59 adjuvant in healthy people of different age groups: a systematic review and meta-analysis Influenza Vaccines for the Future Estimating influenza latency and infectious period durations using viral excretion data Effectiveness of Flu Shot Is 60%-in a Good Year Prevention and control of seasonal influenza with vaccines: recommendations of the Advisory Committee on Immunization Practices-United States, 2019-20 influenza season Intraseason waning of influenza vaccine protection: evidence from the US influenza vaccine effectiveness network This year's flu shot doesn't match what's circulating. Here's what that means Prevalence of immunosuppression among US Adults Influenza in immunosuppressed populations: a review of infection frequency, morbidity, mortality, and vaccine responses The natural history of influenza infection in the severely immunocompromised vs nonimmunocompromised hosts Respiratory tract reinfections by the new human Metapneumovirus in an immunocompromised child The authors have no conflicts of interest to disclose.The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request. Department of Analytics and Operations Management Bucknell University, Lewisburg, PA.