Summary
Seasonal influenza viruses typically cause annual epidemics worldwide infecting 5-15% of the human population1. However, during the first two years of the COVID-19 pandemic, seasonal influenza virus circulation was unprecedentedly low with very few reported infections2. The lack of immune stimulation to influenza viruses during this time, combined with waning antibody titres to previous influenza virus infections, could lead to increased susceptibility to influenza in the coming seasons and to larger and more severe epidemics when infection prevention measures against COVID-19 are relaxed3,4. Here, based on serum samples from 165 adults collected longitudinally before and during the pandemic, we show that the waning of antibody titres against seasonal influenza viruses during the first two years of the pandemic is likely to be negligible. Using historical influenza virus epidemiological data from 2003-2019, we also show that low country-level prevalence of each influenza subtype over one or more years has only small impacts on subsequent epidemic size. These results suggest that the risks posed by seasonal influenza viruses remained largely unchanged during the first two years of the COVID-19 pandemic and that the sizes of future seasonal influenza virus epidemics will likely be similar to those observed before the pandemic.
Main text
The incidence of seasonal influenza has been unusually low since the start of the COVID-19 pandemic in early 20205,6, with cases reported to WHO remaining >80% below historical averages as of January 20222,7. This dramatic reduction is likely due to non-pharmaceutical interventions aimed at reducing transmission and spread of SARS-CoV-28,9, which are also effective in limiting exposure to seasonal influenza viruses. The global lull in influenza virus circulation during the past two years and consequent lack of immune stimulation has led to widespread concerns of increased susceptibility to seasonal influenza viruses in the population due to waning immunity, potentially resulting in larger and more severe epidemics in upcoming seasons9–11. Previous studies of antibody titres to seasonal influenza viruses prior to the COVID-19 pandemic showed that antibody tires against influenza A viruses typically wane to half peak levels 3.5-10 years after infection12–14. However, evidence is lacking on how antibody immunity against seasonal influenza viruses has changed during the near-absence of seasonal influenza in the COVID-19 pandemic and the impact this could have on future influenza epidemic size.
To investigate how the lack of influenza virus circulation since the start of the COVID-19 pandemic has impacted antibody levels against seasonal influenza viruses, we measured antibody titres, based on haemagglutination inhibition (HI), to representative strains of seasonal A/H3N2, A/H1N1pdm09, B/Yamagata, and B/Victoria viruses in longitudinal serum samples collected every summer between 2017 and 2021 from 100 healthy male adults participating in the Amsterdam Cohort Studies on HIV infection and AIDS (ACS)15 (Fig. 1a, Extended Data Fig. 1a, Extended Data Fig. 2, see Methods). Results from virological and syndromic surveillance in the Netherlands showed that influenza A/H3N2, A/H1N1pdm09 and B/Yamagata viruses caused epidemics in the three influenza seasons prior to the onset of the COVID-19 pandemic (Fig. 1b) and that the epidemic activity during this period was consistent with patterns from 2010-2019 (Extended Data Fig. 3). For all seasonal influenza virus types and subtypes, cohort mean HI titres increased after the 2017/2018 influenza epidemic but returned to pre-2017/2018 levels by summer 2019 (Fig. 1a and Extended Data Fig. 2). From 2019 until 2021, mean HI titres remained largely unchanged for all influenza virus (sub)types, including during the COVID-19 pandemic period when there was negligible influenza virus circulation (Fig. 1a and Extended Data Fig. 2). Differentiating the year-on-year individual HI titre distributions by titre rises that are indicative of recent influenza virus infection (≥4-fold increase, ≥2 log2 units), showed that influenza A and B virus infections were most common in individuals with low antibody titres in the year prior to infection (Fig. 1c and Extended Data Fig. 2). Overall, the HI titre distributions of the cohort remained largely unchanged over the study period, including during the COVID-19 pandemic.
a, Individual antibody titres against seasonal influenza viruses based on haemagglutination inhibition (HI) assay from 2017-2021 among 70 healthy male adult participants of the Amsterdam Cohort Studies on HIV infection and AIDS (ACS) cohort for each influenza virus (sub)type (see Methods). Mean antibody tires changes across all individuals are drawn in bold lines with error bars indicating the mean standard error (n=70). b, Seasonal influenza virus epidemic activity 2017-2021 in the Netherlands based on virological and syndromic surveillance data. c, HI titre distributions in the ACS cohort following each winter epidemic period coloured by influenza virus (sub)type. HI titre distributions of individuals who experienced a ≥2 log2 units increase in HI titre (≥4-fold increase in HI titre), indicating likely infection in the previous winter epidemic period, are shown in grey bars. d, HI antibody titre waning rates by influenza virus (sub)type in adults estimated from HI titres from 70 ACS and 65 RECoVERED participants. Error bars correspond to the 50% and 95% credible interval from the Markov Chain Monte Carlo algorithm used to explore the distribution of model parameters. Waning rate of -1.0 corresponds to one two-fold decrease in antibody titre in one year.
Antibody waning during pandemic
Using a mathematical model on the HI titres of ACS participants who were likely not infected during the entire 2017-2021 period, i.e. no ≥2 log2 unit increases in HI titre, we estimated that antibody titres against A/H3N2 viruses waned at -0.20 log2 units per year, 95% credible interval (CI) (−0.24, -0.16); A/H1N1pdm09 viruses at -0.10, 95%CI (−0.12, -0.07); B/Victoria viruses at -0.13, 95%CI (−0.16, -0.10); and B/Yamagata viruses at -0.14, 95%CI (−0.17, -0.11) (Fig. 1d, Extended Data Fig. 2, and Extended Data Tables 1 and 2), in agreement with waning rates previously reported for adults12,14. This indicates that substantial waning of immune protection against seasonal influenza viruses occurs at timescales that are substantially longer than the lull in seasonal influenza virus circulation during the first two years of the COVID-19 pandemic. We also calculated waning rates using HI titres from the same ACS individuals but only for the period after the start of the COVID-19 pandemic, i.e. 2020-2021 (Fig. 1d, Extended Data Fig. 2, and Extended Data Tables 1 and 2). There was also no significant waning of HI titres against any of the viruses during this period.
To extend our observations beyond the ACS cohort, we also measured antibody titres to the same representative viruses in serum samples collected in mid-2020 and mid-2021 from a longitudinal cohort of adult COVID-19 patients who were confirmed not to be vaccinated against seasonal influenza viruses in 2020 (Extended Data Fig. 1b) (Viro-immunological, clinical and psychosocial correlates of disease severity and long-term outcomes of infection in SARS-CoV-2 – a prospective cohort study, referred to as RECoVERED16). In this cohort, we estimated waning rates towards A/H3N2, A/H1N1pdm09, B/Yamagata, and B/Victori to be -0.15, 95%CI (−0.31, 0.01), -0.08, 95%CI (−0.19, 0.03), -0.08, 95%CI (−0.20, 0.04) and - 0.10, 95%CI (−0.22, 0.02) log2 units per year respectively, in good agreement with those derived from the ACS cohort (Fig. 1d, Extended Data Fig. 2, and Extended Data Table 1). Combining data from both cohorts for the 2020-2021 period, the estimated waning rate remained similar to previous estimates for A/H3N2, A/H1N1 and B/Yamagata, and negligible for B/Victoria (Extended Data Table 1). Measurement error was found to be consistent in both datasets at 0.38, 95%CI (0.36, 0.40) and 0.33, 95%CI (0.31, 0.36) log2 units for the full ACS and RECoVERED cohorts respectively (Extended Data Table 1), corresponding to a one-sided probability of a 2-fold error of approximately 6 −10%. Taken together, these results suggest that there have only been negligible changes in antibody titres to seasonal influenza viruses among adults since the start of the COVID-19 pandemic.
The lack of HI antibody titre waning suggests that immunity to seasonal influenza viruses in adults is unlikely to have declined substantially during the first two years of the pandemic. However, previous work showed that waning in children might be different from adults and could have an impact on susceptibility to infection14. While not possible to investigate waning in children in our cohorts, historical lulls in circulation of particular influenza virus (sub)types and their impact on subsequent epidemics have the potential to offer insights into how changes in population immunity, or lack thereof, could impact seasonal influenza epidemics in the post-COVID-19 pandemic period.
Past lulls in influenza virus circulation
Prior to the COVID-19 pandemic, seasonal influenza virus circulation was highly heterogeneous with influenza epidemics typically being dominated by one or two of the four seasonal influenza viruses. This heterogeneity led to frequent periods of 1-3 years where some of the seasonal influenza virus subtypes barely circulated in many countries. These periods are potentially analogous to the first two years of the COVID-19 pandemic as the near-absent circulation of some influenza virus (sub)types might render populations more susceptible, hence leading to larger epidemics.
To investigate the impact of this heterogeneity, we computed the proportion of isolates corresponding to each influenza (sub)type in 718 season-country pairs from 2002/2003 to 2019/2020, for 47 countries in temperate zones (Fig. 2a), based on influenza virus virological surveillance data deposited in the WHO FluNet database7. Low or near-absent circulation of an influenza virus (sub)type within a single season (i.e. a (sub)type accounting for <10% of all influenza virus isolates in a country’s season) occurred regularly during this period accounting for 29%, 27%, and 38% of all country-season pairs for A/H3N2, A/H1N1pdm09 and influenza B viruses respectively (Fig. 2b).
a, Geographic distribution of countries included in the dataset for epidemic composition by subtype, coloured by number of seasons. b, Proportion of viral isolates in each epidemic by virus (sub)type, across all countries and seasons (0 indicates complete absence, 1 indicates complete dominance). c, Geographic distribution of countries included in the dataset for epidemic size by (sub)type. d, Relative epidemic sizes by virus (sub)type, across all countries and seasons (the lower the number, the smaller the epidemic). e, Epidemic dominance as a function of years since previous dominance, by (sub)type. Error bars correspond to 95% confidence interval from an exact two-tailed binomial test for proportions. f, Relative size of a (sub)type-specific epidemic as a function of number of years since previous dominance of that (sub)type in that country, coloured by season. Each point corresponds to a country. g, Relative size of a subtype’s epidemic as a function of its size in the previous season and the sum of the two previous seasons’ sizes. Each point corresponds to a country-season, coloured by the season. h, Posterior distributions of parameter estimates in a Bayesian hierarchical model for epidemic size, with one year since previous dominance (circles), previous epidemic size (diamonds), or sum of previous two epidemics’ size (squares) as predictors, either with or without season effects. Estimates above the dotted line represent the season effects. Points, thick and thin lines correspond to the posterior mean, 50% CI, and 95% CIs, respectively.
While virological surveillance data yield insights into the frequency of circulation of influenza viruses, it can be biased for estimating epidemic size due to year-on-year changes in virus sampling rates. To complement the virologically confirmed dataset used above, we estimated type and subtype-specific epidemic size based on influenza-like illness (ILI) data from the WHO FluID17 database for 20 countries in Europe and the Middle East in the period from 2010 to 2020 (Fig. 2c). To compute epidemic sizes for each influenza (sub)type in each season, we multiplied a country’s ILI burden by the proportion of isolates attributed to each type and subtype in each season. We then divided this number by the total ILI burden across all ten seasons and multiplied this number by the total number of seasons to calculate relative epidemic sizes. In these estimates, a relative size of one corresponds to the mean number of influenza virus infections in a single season for a given country, irrespective of type and subtype. Epidemic sizes lower than 0.1, indicating very small or absent epidemics, were observed for 28%, 23%, and 37% of country-seasons for A/H3N2, A/H1N1pdm09 and influenza B viruses, respectively (Fig. 2d).
To investigate the effect of periods of low influenza virus circulation on epidemic (sub)type composition, we calculated the probability of an influenza virus (sub)type’s dominance as a function of years since previous dominance, where we defined dominance as a (sub)type accounting for >30% of a season’s isolates (Fig. 2e). The probability of a (sub)type’s dominance increased with greater number of years since previous dominance. However, this analysis also implies that there were periods of up to three years where an influenza (sub)type did not dominate in the past, indicating that periods of low to absent circulation of particular seasonal influenza viruses were also common before the COVID-19 pandemic. Mean epidemic sizes for each influenza virus (sub)type increased with time since dominance (Fig. 2f). However, these increases were strongly related to probability of (sub)type dominance (Fig. 2e) and epidemic sizes varied substantially since last dominance, suggesting that the overall influence of time since dominance on epidemic size is relatively small.
Epidemic sizes of each (sub)type have a negative relationship with incidence of that specific (sub)type in the preceding year with large successive epidemics being rare (Fig. 2g left column). However, this negative relationship largely disappears when taking into account the cumulative incidence of each (sub)type two years into the past (Fig. 2g, right column), with epidemics of high and low incidence both being likely to occur if preceded by years of low-to-mid incidence.
For each number of years since dominance, there is a striking degree of clustering of epidemic size across countries by season (Fig. 2f). For example, in the 2013/2014 and 2016/2017 seasons, where A/H3N2 dominated in most countries two years prior, the relative incidence in 2016/2017 appeared consistently higher than in 2013/2014. Notably, in 10 of the 20 countries included in our dataset, the first A/H3N2-dominant season (2011/2012) after the 2009 A/H1N1pdm09 pandemic did not belong to the three largest A/H3N2 epidemics in the influenza seasons from 2010-2011 until 2019-2020, despite three years of near absent circulation.
Season effects dominate epidemic size
Given the degree of clustering of epidemic size by season, we hypothesized that the size of influenza virus (sub)type-specific epidemics could be explained by a combination of season-specific effects and effects related to the presence or absence of that virus (sub)type in the years preceding an epidemic. Season-specific effects, shared between countries in a single season, may include a variety of viral, host, environmental, or epidemiological variables, such as antigenic novelty, climate, heterosubtypic competition, or the flux of viral seeding18–21
To investigate this hypothesis, we constructed a Bayesian hierarchical model that uses these effects as predictors of the (sub)type-specific size of seasonal influenza epidemics. Specifically, we considered models that individually include the number of years since previous dominance of that (sub)type, the size of that (sub)type’s epidemic in the previous year, or the sum of that (sub)type’s incidence in the previous two years as predictors. In these models, the season effects correspond to the predicted ‘base size’ of a (sub)type’s epidemic, given that the previous dominance was in the previous year, given that there was no circulation in the previous year, or given that there was no circulation in the previous two years, respectively. These season effects are modulated by the effects of prior circulation to yield an epidemic’s predicted size. Years since dominance, size in the previous year and the sum of previous two seasons’ sizes individually had non-trivial effects on epidemic size (Fig. 2h). However, between-season differences with regard to season effects were consistently of substantially greater magnitude than any of the predictors related to prior incidence across all model formulations, suggesting that season-specific factors unrelated to the absence or presence of viral circulation in the previous year(s) dominate epidemic size. Previous epidemic size appeared to have a moderate effect on epidemic size. This effect substantially decreased when looking at the sum of the two previous epidemic sizes (Fig. 2h).
Because epidemic sizes clustered strongly by season, which might obscure the effect of prior incidence in a model with season effects, we also considered models without season effects. Here, estimates on the impact of absence or presence of circulation in prior years on size were higher, but the differences between seasons with regard to season effects (‘base sizes’) remained far greater (Fig. 2h). For example, even when using parameters estimated from a model without season effects, the model predicts that the size of an A/H3N2 epidemic with the mean estimated season effect and previous dominance three years prior is smaller than an epidemic with the largest estimated season effect and A/H3N2 domination in the previous season. Models that included season effects exhibited much better predictive performance than models without season effects (Extended Data Fig. 4a). Additionally, models that included prior incidence of the opposite subtype had substantial effects of opposite sign, implying that the estimated effects of prior incidence might reflect a combination of prior incidence and effects of heterosubtypic competition (Extended Data Fig. 4b). These results suggest that inherent season-specific effects have more substantial effects on epidemic size than (sub)type-specific patterns of prior circulation.
Taken together, the lack of changes observed in the pattern of measured antibody titres against seasonal influenza viruses and nearly two decades of epidemiological data suggest that the near-absence of seasonal influenza virus circulation during the first two years of the COVID-19 pandemic is unlikely to result in substantially larger influenza epidemics in the years to come. The size of future influenza epidemics is likely to fall within the size distribution of epidemics in the years before the COVID-19 pandemic.
Methods
Viruses
Based on phylogenetics data, four seasonal influenza viruses, A/Netherlands/04189/2017 (H3N2), A/Netherlands/10218/2018 (H1N1), B/Netherlands/04136/2017 (Yamagata), and B/Netherlands/00302/2018 (Victoria), were selected as representatives of the four main types/subtypes of seasonal influenza viruses that circulated prior to the SARS-CoV-2 pandemic. To select these strains, we downloaded high-quality (<5% ambiguous nucleotides, >95% full length) seasonal influenza virus haemagglutinin sequences (A/H3N2, N=1,396; A/H1N1pdm09, N=1,283; B/Yamagata, N=1,129; and B/Victoria, N=1,408) collected between 2016 and October 2021 from GISAID (www.gisaid.org) and reconstructed maximum-likelihood phylogenetic trees for each influenza virus subtype using the general time reversible substitution model with IQ-TREE22. These trees were used to assess the representativeness of viruses from the Netherlands in the early portion of the study period and the selected viruses were all representative of viruses that caused epidemics in the Netherlands during the 2017/18 winter. All four viruses were propagated on Madin-Darby Canine Kidney (MDCK) cells in infection medium, which consisted of MEM-Eagle Medium /EBSS (Lonza, Geleen, The Netherlands) supplemented with MEM Non-Essential Amino Acids (Gibco, ThermoFischer Scientific, Amsterdam, The Netherlands), penicillin (100 U/mL), streptomycin (100 g/mL), L-Glutamine (Lonza), HEPES (Lonza), and TPCKtrypsin (Sigma-Aldrich/Merck, Darmstadt, Germany). They were harvested after 72 hours of incubation at either 37°C (H3N2 and H1N1) or 33°C (Yamagata and Victoria) and checked by Sanger sequencing.
Longitudinal serum samples
A total of 630 serum samples from 165 healthy adults were collected in the Netherlands, longitudinally before and during the pandemic in two separate cohorts: 1. Amsterdam Cohort Studies on HIV infection and AIDS15 (ACS) (100 participants and a total of 500 samples) and 2. the Viro-immunological, clinical and psychosocial correlates of disease severity and long-term outcomes of infection in SARS-CoV-2 – a prospective cohort study16 (RECoVERED) (65 participants and a total of 130 samples).
The initial aim of the Amsterdam Cohort Studies was to investigate the prevalence, incidence, and risk factors of HIV-1 infection. The study population consists of men who have sex with men and live mainly around the city of Amsterdam, the Netherlands. Enrolled men were both HIV-1 seronegative and HIV-1 seropositive. Participants from the ACS cohort included in our study were all HIV-1 seronegative men ranging from 22 to 70 years old at the time of first collected sample. Briefly, five samples were collected per participant, i.e. 1. mid-2017, 2. mid-2018, 3. mid-2019, 4. mid-2020, 5. mid-2021. Using only samples collected in the summer period potentially helps to overcome issues that could arise from the transient antibody boosts due to both infection and vaccination23. The Amsterdam Cohort Studies on HIV infection and AIDS was approved by the Medical Ethics Committee of the Amsterdam University Medical Centre of the University of Amsterdam, the Netherlands (MEC 07/182). Participation in ACS is voluntary and without incentive. Written informed consent of each participant was obtained at enrolment.
Data derived from the ACS samples was complemented by the data acquired from participants in the SARS-CoV-2 cohort RECoVERED. The aim of the RECoVERED cohort study is to describe the immunological, clinical and psychosocial sequelae of a SARS-CoV-2 infection. Individuals aged 16 to 85 years with laboratory-confirmed SARS-CoV-2 infection were enrolled on May 2020 till the end of June 2021 in the municipal region of Amsterdam, the Netherlands. The RECoVERED study was approved by the medical ethical review board of the Amsterdam University Medical Centre (NL73759.018.20). All participants provided written informed consent. From RECoVERED, we selected a total of 65 individuals, both male and female adults ranging from 20 to 77 years old at the time of the first collected sample, all of which had a confirmed SARS-CoV-2 infection but were otherwise healthy and unvaccinated for influenza in 2020. For these 65 individuals, samples were collected in the summer period of 2020 and 2021 only (two total for each participant). Samples from the RECoVERED cohort add diversity in age and gender to the ACS samples set but lacked pre-pandemic samples.
Haemagglutination Inhibition (HAI) assay
All serum samples were receptor destroying enzyme (RDE)-treated, as described elsewhere24. Briefly, 100mL of serum samples from ACS individuals 1-30 were combined with 200mL of RDE; for serum samples from ACS individuals 31-100 and all 65 RECoVERED subjects, 100mL of serum were combined with 300mL of RDE. This difference in protocol was per the instructions of the providers of the respective batches of RDE. Because of this protocol difference, the results of ACS participants 1-30 and 31-100 are shown separately. All samples were then incubated at 37°C for 18 to 20 hours. The RDE reaction was then halted by heating the treated samples at 56°C for 30 to 60 minutes.
The haemagglutination inhibition activity of all serum samples was tested in an haemagglutination inhibition assay as described elsewhere24,25 using two replicates per sample for A/H1N1, B/Yamagata, and B/Victoria, and one single measurement for A/H3N2. Briefly, the haemagglutination titre of each of the four viruses was determined by doing a two-fold serial dilution of 50mL of each virus stock and adding 50mL of PBS and 25mL of 1% turkey red blood cells (tRBCs) to each well, followed by one hour incubation at 4°C and the reading of the haemagglutination patterns. The virus stocks were then diluted to a concentration of 4 haemagglutination units (HAU). The diluted viruses were then incubated with 50mL of two-fold serially diluted serum, in a total volume of 75mL for 30 minutes at 37°C. The initial dilution used for the serial dilution of the serum was 1:20 of the RDE treated serum. After the incubation step, 25mL of 1% turkey red blood cells were added to the serum-virus mix and incubated at 4°C for one hour. The haemagglutination inhibition patterns were then read out and used for the calculation of antibody titres. Due to the known inefficient agglutination of tRBCs by recent A/H3N2 viruses, we used glycan remodelled turkey red blood cells expressing appropriate receptors for recent A/H3N2 viruses26 for the implementation of the assay for the A/H3N2 virus stock.
Data pre-processing
Two approaches were used in selecting data from which to determine antibody waning rates. Firstly, we used samples from the RECoVERED cohort for the years 2020 and 2021, where all participants are confirmed to have not received an influenza vaccination between the two sample collections and no natural influenza infection can be safely assumed given the near absence of influenza in the Netherlands during this period. Second, we used the ACS data for the 5 years from 2017 to 2021. Individuals who experienced a 4 or greater fold increase in titre between consecutive visits for a particular strain had their data for the strain removed in order to remove the obscuring effects of vaccination and infection. The advantage of the former approach is the certainty regarding infection and vaccination status, the latter, however, allows a longer period of time over which to observe potentially subtle antibody waning dynamics.
Antibody waning model
True antibody titre log2 HI, 
, as opposed to that measured by haemagglutination inhibition assay, Ti, is a continuous variable which we assume, for every individual, i, decays with time, t, as
Where ci are individual specific initial titres and α is the shared waning rate.
If serum dilutions could be performed in arbitrarily small increments, we assume the point at which haemagglutination would be observed to cease, Tobs, to be distributed normally about the true value:
where ϵ shall be referred to as the “measurement error”. Instead, with discrete dilutions in increments of 1, the probability of measuring T ∈ {0, 1, 2…8} is the probability that Tobs falls between T and T − 1. Thus, the measurement probability is given by:
where Φ(x, μ, σ) is the cumulative distribution function of the normal distribution.
Our data for each individual, i, consists of a series of titre measurements, 
, at corresponding timepoints 
, where r ∈ {1,2} indicates replicate measurements. To infer the probability of the unknown parameters ϵ and a given the data, it is necessary to augment the data by introducing individual intercepts. For one replicate from one individual, the likelihood of unknown parameters a, ϵ, and ci then becomes:
where 
are the true values of titre, given the unknown parameters and Π is the prior joint distribution of the parameters. The total log likelihood is thus the sum over all individuals and replicates:
A Markov Chain Monte Carlo (MCMC) algorithm was used to explore the distribution of model parameters (waning and measurement) and augmented data (individual intercepts). In each iteration, model parameters were updated using a Metropolis-Hastings (MH) algorithm and 10% of participants were randomly selected for updating of their augmented data, also via MH. This model was run on 4 independent chains, each consisting of 20,000 iterations discarding the first 5,000 as burn in. Non-informative priors were used and convergence was assessed by inspection of the trace plots and Rhat from Stan v2.21.0. Analyses were conducted using R v4.0.3, with code available in the Github repository.
Epidemic composition data
We downloaded records of virological surveillance data from the WHO FluNet7 database for all countries in the temperate Northern and Southern Hemisphere from 2002 until 2020, or a shorter period for a limited subset of countries. For each country, we retained the longest sequence of consecutive seasons in which at least 20 specimens were influenza-positive. In each season, defined as the period from the period from week 40 until week 20 for the Northern Hemisphere and the entire year for the Southern Hemisphere, we computed the proportion of all positive tests that was attributable to each of A/H3N2, pandemic A/H1N1pdm09 (from the 2009 pandemic onwards), and influenza B viruses. We did not break down influenza B viruses by lineage because in many countries influenza B viruses were not further characterized. In many seasons, only a proportion of all influenza A virus positive tests were subtyped; in those cases, we approximated the total proportion of each subtype by assuming that the subtype of the non-subtyped influenza A virus specimens were distributed according to the relative proportions of subtyped influenza A viruses. This resulted in a dataset of 679 season-country records over a period of 18 seasons, for 46 countries. We assigned a binary variable to each subtype/type in each season for dominance; we defined dominance as a subtype/type accounting for at least 30% of all isolates in a country in a season. Hence, in principle, all three subtypes/types considered can be simultaneously dominant in a single season. To avoid including effects of the COVID-19 pandemic on influenza dynamics, we truncated the 2019-2020 season at the 15th of February 2020, and to avoid including the effect of the 2009 A/H1N1pdm09 pandemic, we truncated the 2008-2009 influenza season at the 1st of April 2009.
Epidemic size data
To estimate epidemic sizes for each subtype/type, we extracted weekly records of influenza-like illness from the WHO FluID17 database. We limited this dataset to countries for which ILI records were available for all seasons from 2010-2011 until 2019-2020, and for which virological surveillance data was available, as described above. Additionally, we required ILI curves to follow the expected shape of an influenza epidemic curve, i.e. peaking in winter and only sporadic isolation outside this period, and without periods of missing data. This yielded a set of 22 countries, each with 10 seasons worth of ILI data. We approximated the epidemic size of each subtype in each country per season by multiplying total ILI incidence in that country’s season by the proportion of all isolates from that country in that season attributable to that subtype in the virological surveillance data. We computed the relative size of each subtype’s epidemic in each country by computing the proportion of all ILI in the period from the 2010-2011 until 2019-2020 seasons that was attributable to that subtype and multiplying this number by the total number of seasons. Because this metric averages across all ten seasons, it gives an estimate of the size of an epidemic in a country, relative to all other seasons in that country, for each subtype. We accounted for the COVID-19 pandemic as described above.
Statistical modelling
To estimate the effect of presence or absence of influenza virus circulation in the previous year(s) on epidemic size, we performed Bayesian hierarchical linear regression. The first model has relative size as outcome and time since dominance as calculated using the virological surveillance data as predictor (N=188, 198, 180 for A/H3N2, A/H1N1pdm09 and B):
where yi is an epidemic’s relative size for a certain (sub)type in country-season pair i, as[i] is the season effect for that (sub)type corresponding to that season, β is the coefficient for number of years since dominance of the (sub)type, and sy is the error standard deviation. xi represents the number of years since previous dominance of the (sub)type in country-season pair i minus one, such that as[i] represents the predicted size if the previous dominance was in the previous year. Season effects are shared between countries in a single season.
We assumed that the season effects as, constrained to positive values, are distributed according to a common mean µa and common standard deviation σa.
We put weakly informative priors on the standard deviations, the mean season effect and the main effect.
We also ran a model with relative size as outcome and relative size in the previous year as the predictor (N=180 for each (sub)type) or the sum of the two previous years size (N=160 for each (sub)type). We use the same model specification and priors as for the size-years since dominance model, but we replace the predictor with the relative size in the previous season, or with the sum of relative size in the two previous seasons. For A/H3N2 and A/H1N1pdm09, we also ran these models with time since dominance, previous season epidemic size and sum of two previous seasons epidemics sizes of the other subtype as predictor (N=198, 188 for A/H3N2, A/H1N1pdm09 for time since dominance, N= 180, 160 for each subtype for previous year, previous two years’ sum, respectively). For all these above models, we also ran the same models without season effects, i.e. with a single value for the intercept, for which the prior is equal to the prior of the mean season effect in the model with season effects.
These models were each run for each subtype individually, for 3000 iterations, discarding the first 1000 as burn-in, with four independent chains. All models were fit using MCMC in Stan v2.21.0, with convergence assessed by inspection of Rhat (< 1.05), effective sample size (> 200) and the trace plots. We compared models with and without season effects using leave-one-out cross-validation27.
Statistics
No statistical method was used to predetermine sample size. Data were not randomized nor analysed in a double-blinded manner. The haemagglutination inhibition activity of all serum samples was tested in an haemagglutination inhibition assay using two replicates per sample for A/H1N1, B/Yamagata, and B/Victoria, and one single measurement for A/H3N2. Specific details regarding amount of data are indicated in the figure captions. For parameter estimates 95% credible intervals were considered as the significant bounds and were calculated from the 2.5th and 97.5th percentiles of the MCMC traces. For Fig. 2e, error bars correspond to 95% confidence interval from an exact two-tailed binomial test for proportions.
Reporting summary
Further information on research design is available in the Nature Research Reporting Summary linked to this paper.
Data availability
Accession codes for GISAID data used in this paper are provided as supplementary information files. Raw de-identified hemagglutination inhibition data and raw surveillance data downloaded from WHO FluNet and FluID can be found in the project GiHub repository (https://github.com/AMC-LAEB/waning-immunity-to-flu). Biological materials are available for study via the Amsterdam Cohort Studies on HIV infection and AIDS (ACS) and the Viro-immunological, clinical and psychosocial correlates of disease severity and long-term outcomes of infection in SARS-CoV-2 – a prospective cohort study (RECoVERED).
Code availability
Custom scripts used for data analysis and modelling are available at GitHub https://github.com/AMC-LAEB/waning-immunity-to-flu.
Author contributions
Z.C.F.G., S.P.d.J., J.G., A.X.H., D.E., M.D.d.J, and C.A.R. designed the research; Z.C.F.G, S.v.L., M.v.H., K.d.H., and L.v.G executed the experimental work; Z.C.F.G. and S.v.L. generated the antibody titre data; E.W., G.J.d.B, H.D.G.v.W., A.M., L.v.d.H., M.P., N.K., and M.D.d.J. collected the clinical samples; R.P.d.V. and G.J.B. made and provided the glycan remodelled turkey red blood cells; S.P.d.J. and J.G. implemented the modelling work and performed the data analysis; Z.C.F.G., S.P.d.J., J.G., A.X.H., B.E.N., and C.A.R. wrote the first draft of the paper. All authors contributed to the critical revision of the paper.
Competing interests
No competing interests.
Acknowledgements
A.X.H., Z.C.F.G. and C.A.R. were supported by ERC NaviFlu (No. 818353). J.G. and C.A.R. were supported by NIH R01 (5R01AI132362-04). C.A.R. was also supported by an NWO Vici Award (09150182010027). R.P.dV. was supported by ERC starting grant 802780 and a Beijerinck Premium of the Royal Dutch Academy of Sciences. GJB was supported by the Netherlands Organization for Scientific Research (NWO TOPPUNT 718.015.003) and by an ERC advanced grant (101020769). The RECoVERED cohort is supported by NWO ZonMw (No. 10150062010002). The Amsterdam Cohort Studies on HIV infection and AIDS, a collaboration between the Public Health Service Amsterdam, the Amsterdam UMC of the University of Amsterdam, Medical Center Jan van Goyen and the HIV Focus Center of the DC-Clinics, are part of the Netherlands HIV Monitoring Foundation and financially supported by the Center for Infectious Disease Control of the Netherlands National Institute for Public Health and the Environment. We gratefully acknowledge the authors and originating and submitting laboratories (supplementary information) for the reference sequences retrieved from GISAID’s EpiFlu Database used in this study. The authors thank Mr. Reinier van der Palen of the department of Chemical Biology and Drug Discovery, Utrecht University for practical assistance with turkey erythrocyte glycan remodeling.
References
Methods references
Subject Area
Reviews and Context
1
Comment
0
TRIP Peer Reviews
0
Community Reviews
0
Automated Services
2
Blogs/Media
Author Videos