key: cord-0948657-ravvlg8l
authors: Domenech de Cellès, Matthieu; Casalegno, Jean-Sebastien; Lina, Bruno; Opatowski, Lulla
title: Estimating the impact of influenza on the epidemiological dynamics of SARS-CoV-2
date: 2021-12-03
journal: PeerJ
DOI: 10.7717/peerj.12566
sha: 55ad243f2f0fa9af431beb4804577be56b56eb8c
doc_id: 948657
cord_uid: ravvlg8l

As in past pandemics, co-circulating pathogens may play a role in the epidemiology of coronavirus disease 2019 (COVID-19), caused by the novel severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2). In particular, experimental evidence indicates that influenza infection can up-regulate the expression of ACE2—the receptor of SARS-CoV-2 in human cells—and facilitate SARS-CoV-2 infection. Here we hypothesized that influenza impacted the epidemiology of SARS-CoV-2 during the early 2020 epidemic of COVID-19 in Europe. To test this hypothesis, we developed a population-based model of SARS-CoV-2 transmission and of COVID-19 mortality, which simultaneously incorporated the impact of non-pharmaceutical control measures and of influenza on the epidemiological dynamics of SARS-CoV-2. Using statistical inference methods based on iterated filtering, we confronted this model with mortality incidence data in four European countries (Belgium, Italy, Norway, and Spain) to systematically test a range of assumptions about the impact of influenza. We found consistent evidence for a 1.8–3.4-fold (uncertainty range across countries: 1.1 to 5.0) average population-level increase in SARS-CoV-2 transmission associated with influenza during the period of co-circulation. These estimates remained robust to a variety of alternative assumptions regarding the epidemiological traits of SARS-CoV-2 and the modeled impact of control measures. Although further confirmatory evidence is required, our results suggest that influenza could facilitate the spread and hamper effective control of SARS-CoV-2. More generally, they highlight the possible role of co-circulating pathogens in the epidemiology of COVID-19.

The current pandemic of coronavirus disease 2019 , caused by the novel severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), has led to global alarm. Following the first case reports in December 2019 in Wuhan, China (Zhu et al., 2020) , SARS-CoV-2 rapidly spread across the globe and has resulted in approximately 171 million cases and 3.6 million deaths worldwide, as of June 1, 2021 (Johns Hopkins University Center for Systems Science and Engineering (JHU CCSE), 2020) . Because of the initial lack of prophylactic or therapeutic treatments, the pandemic caused the implementation of unprecedented control measures, which culminated in the lockdown of several billion people in over 100 countries during April-May 2020 (Hale et al., 2020) . Although a number of fixed (e.g., greater age, male sex) and chronic (e.g., hypertension, diabetes) risk factors of mortality have now been identified (Williamson et al., 2020) , the time-varying drivers of COVID-19 epidemiology remain poorly understood. Experience gained from past pandemics has highlighted the potentially large contribution of co-circulating pathogens to the burden of an emerging disease (Morens, Taubenberger & Fauci, 2008) . Despite the relevance for epidemic forecasting and for designing control strategies, the impact of co-circulating pathogens on SARS-CoV-2 epidemiology has remained largely unexplored (Thindwa et al., 2020) .

Respiratory viruses-including SARS-CoV-2 and other coronaviruses, rhinoviruses, influenza viruses, etc.-form a large class of viruses that cause seasonal infections of the respiratory tract in humans. Mounting evidence indicates that their epidemiologies are not independent, as a result of interaction mechanisms that may operate at different scales and that can be classified as either facilitatory or antagonistic (DaPalma et al., 2010; Opatowski, Baguelin & Eggo, 2018) . The interaction between the respiratory syncytial virus (RSV) and influenza may provide an example of antagonism. Indeed, experimental evidence in ferrets has shown that influenza viruses induce an antiviral state that transiently limits secondary infection with RSV (Chan et al., 2018) , an effect postulated to explain the delayed epidemic of RSV during the 2009 influenza pandemic (Casalegno et al., 2010; Mak et al., 2012) . Although such antagonistic interactions appear, to date, to be the most common among respiratory viruses (Opatowski, Baguelin & Eggo, 2018) , experimental evidence indicates that co-infections may also have a facilitatory effect, for example by increasing viral growth (Goto et al., 2016) . Increased transmission of influenza during co-infection with other respiratory viruses was also proposed to explain the multiple waves during the 1918 influenza pandemic (Merler et al., 2008) . Interestingly, according to recent evidence a viral respiratory infection (in particular with influenza viruses) can up-regulate the expression of ACE2-the cognate receptor of SARS-CoV-2 in human cells-in the respiratory epithelium (Smith et al., 2020; Ziegler et al., 2020) . In addition, this up-regulation was demonstrated experimentally to increase infectivity of SARS-CoV-2 in mice co-infected with influenza A virus (Bai et al., 2021) . This suggests that respiratory viruses could affect the epidemiology of SARS-CoV-2. Here, we hypothesized that influenza-which peaked in February 2020 and therefore co-circulated during the early spread of COVID-19 in Europe (Fig. 1B) -interacted with SARS-CoV-2.

To test this hypothesis, we developed a semi-mechanistic, population-based model of SARS-CoV-2 transmission and of COVID-19 mortality. Using likelihood-based statistical inference methods, we fitted this model to mortality incidence data in four European countries to test a range of assumptions about the impact of influenza on the transmission (Flaxman et al., 2020) . (B) Time plot of influenza incidence, calculated as the product of the incidence of influenza-like illnesses and of the fraction of samples positive to any influenza virus (see also Fig. S1 for a time plot of the latter two variables). The vertical dashed lines delimitate the period of overlap between SARS-CoV-2 and influenza, defined as the period between the assumed start date of SARS-CoV-2 community transmission and 6 weeks after the epidemic peak of influenza . In each country, the time series displayed were incorporated as covariates, which modulated the transmission rate of SARS-CoV-2 in our model (see Methods). In B, the y-axis values differ for each panel.

Full-size DOI: 10.7717/peerj.12566/ fig-1 dynamics of SARS-CoV-2. We find that influenza may have transiently increased the transmission of SARS-CoV-2 during the first wave of COVID-19 in Europe. Although further confirmatory evidence is required, our results suggest that influenza could facilitate the spread of SARS-CoV-2 and, more generally, emphasize the potential role of cocirculating pathogens in the epidemiology of COVID-19.

Country-level time series of the stringency index were available from the Oxford COVID-19 Government Response Tracker, developed at the University of Oxford and described elsewhere (Hale et al., 2020) . Briefly, the stringency index provides an aggregate measure of the number and of the strictness of non-pharmaceutical control measures implemented by governments in response to the COVID-19 epidemic. The stringency index is defined as the average of nine normalized ordinal variables, which quantify the strength (e.g., recommended or required) and the scope (e.g., targeted or general) of closure and containment measures (8 variables) and of health measures (1 variable). The resulting index allows quantifying the strength of control measures in a systematic way, on a scale ranging from 0 (no interventions) to 100 (maximum number and maximal intensity of control measures). Of note, however, the stringency index does not quantify the impact of control measures, which likely varied across countries (Flaxman et al., 2020) . In formulating our model, we therefore modeled the relationship between the stringency index and the relative reduction in SARS-CoV-2 transmission using a non-decreasing function, whose parameters represented the impact of control measures and were estimated from the data.

Virological data on the weekly numbers of samples tested and of samples positive to any influenza virus were available from the FluNet database, compiled by the WHO (Fig. S1A ). Parallel syndromic data on the weekly incidence rate of influenza-like illnesses (ILI) were available from the FluID database, also compiled by the WHO (Fig. S1B ). These data were deemed high-quality and used in a previous study on influenza forecasting in the countries considered here (Kramer & Shaman, 2019) . The weekly incidence rate of influenza was then calculated as the product of ILI incidence and of the fraction of samples positive to any influenza virus (Fig. 1B) . Because the magnitude of influenza incidence thus calculated varied markedly across countries (e.g., as a result of different surveillance systems and case definitions), we rescaled each time series by its average during the period of co-circulation of influenza and SARS-CoV-2 (Fig. 1B) . The resulting time series was therefore dimensionless and equalled 1 when influenza incidence equalled its average value during that period.

Data on the daily number of deaths caused by SARS-CoV-2 (counted by date of death) were available from national public health public institutes, in Belgium (Belgian institute for health, sciensano)(2020), in Italy (Dipartimento della Protezione Civile, 2020), and in Spain (Instituto de Salud Carlos III, official data with historical corrections compiled by the media (DATADISTA, 2020). In Norway, the data were available from the worldwide database compiled by the European Center for Disease Control and Prevention (European Centre for Disease Prevention and Control, 2020) . Following a previous study (Flaxman et al., 2020) , and to avoid a possible bias caused by the dominance of deaths due to non-locally acquired infections early in the epidemic, we included observed deaths from the date after which the cumulative observed death count exceeded 10. Data points before that date were treated as missing and were assigned a conditional log-likelihood of 0, such that they did not contribute to the overall log-likelihood. The data were not further pre-processed, except in Italy, where a negative death count was reported on 24 June 2020 and was treated as missing and also assigned a log-likelihood of 0.

We formulated a variant of the standard Susceptible-Exposed-Infected-Recovered transmission model (Keeling & Rohani, 2008) , using the method of stages to allow for realistic distributions of the latent, infectious, and onset-to-death periods (Lloyd, 2001; Wearing, Rohani & Keeling, 2005) . Specifically, we assumed that the latent and infectious periods were Erlang-distributed with shape parameter 2 and mean 1/σ = 4 days and 1/γ = 5 days, respectively (Li et al., 2020) . The resulting generation time T g (i.e., the time from infection of a primary case to transmission to a secondary case) had a mean of 6.5 days and a coefficient of variation of 0.58 (see Fig. S2 for the full distribution and the details of the calculation), to Svensson (2007) and Camacho et al. (2011) consistent with empirical observations and with the values fixed in a previous modeling study (Flaxman et al., 2020; Bi et al., 2020) . To model the impact of the gradual implementation of nonpharmaceutical control measures (e.g., border closure, school closure, lockdown), we mapped the stringency index (denoted by si(t )) to the time-varying relative reduction in transmission of SARS-CoV-2 (denoted by r β (t )). Specifically, we used the following simple linear scaling function, with saturation:

Here the parameter b quantifies how steeply the transmission rate of SARS-CoV-2 decreases as the stringency index increases. Hence, this parameter can be interpreted as a measure of the impact of non-pharmaceutical control measures on SARS-CoV-2 transmission. The deterministic variant of the model was represented by the following set of differential equations:

The force of infection (that is, the per capita rate at which susceptible individuals contract infection (Keeling & Rohani, 2008) ), λ(t ), was modeled as:

where R 0 represents the basic reproduction number of SARS-CoV-2, β 0 = R 0 γ the basic transmission rate, R e (t ) the time-varying effective reproduction number, N the population size (assumed constant during the study period), and F (t ) the renormalized time series of influenza incidence, incorporated as a covariate into the model (Fig. 1) . With this formulation, the parameter β F quantifies the impact of influenza on SARS-CoV-2 transmission: β F > 0 if influenza increases transmission, β F < 0 if influenza decreases transmission, and β F = 0 if influenza has no impact on transmission (null hypothesis). More specifically, the average incidence of influenza during the period of co-circulation with SARS-CoV-2 corresponds to F (t ) = 1, such that 1+β F represents the average relative variation of SARS-CoV-2 transmission associated with influenza. In writing the equations, we implicitly assume that the impact of influenza on SARS-CoV-2 transmission, if any, is short-lived and does not extend long after influenza infection. Finally, we incorporated an observation model that related the dynamics of SARS-CoV-2 infection to that of COVID-19 mortality, taking into account the fact that only a fraction of infections results in death and that, among those, death occurs some time after symptom onset (Flaxman et al., 2020; Verity et al., 2020; Khalili et al., 2020) . We assumed an average duration of pre-symptomatic of 2.5 days, resulting in an average incubation period of 6.5 days, in broad agreement with previous empirical studies ( Khalili et al., 2020; Tindale et al., 2020) . Hence, individuals in the first infected state (I 1 ) were considered pre-symptomatic, and the onset of symptoms was assumed to coincide with the transition from I 1 to I 2 . The onset-to-death time was then assumed to be Erlang distributed with shape parameter 5 and mean 1/κ = 17.8 days (coefficient of variation of 0.45), the value estimated in a previous epidemiological study (Verity et al., 2020) . In a sensitivity analysis, we also tested a mean onset-to-death time of 1/κ = 13 days, the lower bound estimated in a meta-analysis (Khalili et al., 2020) . According to previous studies in European countries, the infection fatality ratio (IFR) typically ranged from 0.5% to 1% early in the epidemic (Salje et al., 2020; O'Driscoll et al., 2021; . We fixed the IFR to µ = 0.01 in the base model, but we considered an alternative value of 0.005 in a sensitivity analysis. Given those assumptions, the observation model was modeled by the following set of ordinary differential equations:

Here D M is the simulated number of daily deaths, modeled as an accumulator variable and reset to 0 at the end of each day. The observed number of daily deaths, D O , was modeled using a negative binomial distribution with mean D M and over-dispersion

, a standard distribution used in previous modeling studies (Flaxman et al., 2020; King et al., 2015) .

As in Flaxman et al. (2020) , simulations were started 30 days before the date from which the cumulative observed death count first exceeded 10. At that date, we assumed that E 1 (0) individuals had been exposed to SARS-CoV-2; other individuals were assumed susceptible to infection (i.e., S(0) = N − E 1 (0)), and all other compartments were initialized to 0.

The stochastic variant of the model was implemented as a continuous-time Markov process approximated via a multinomial modification of the τ -leap algorithm (He, Ionides & King, 2010) , with a fixed time step of t = 10 −1 day. To model the effect of superspreading events-a key feature of SARS-CoV-2 transmission dynamics (Althouse et al., 2020)-, extra-demographic stochasticity was added to the transmission rate β(t ). Specifically, as proposed by Kain et al. (2020) , at every time step we drew a value of β 0 from a Gamma white noise distribution:

with mean µand variance µσ 2 . Here k represents the dispersion parameter of the Negativebinomial distribution for the individual reproduction number (with mean R 0 and variance

as estimated in previous studies (Lloyd-Smith et al., 2005; Endo et al., 2020) . As in Kain et al. (2020) and in keeping with empirical estimates from contact tracing studies of SARS-CoV-2 (Endo et al., 2020) , we fixed k = 0.16.

Following the method presented in Flaxman et al. (2020) , we estimated unknown model parameters using observed COVID-19 mortality data alone. Indeed, because of the initially limited testing capacity (typically reserved to severe cases or high-risk groups), mortality data were arguably more reliable than case data early in the epidemic in most countries (Flaxman et al., 2020) . By incorporating known epidemiological parameters (key among those the onset-to-death time, the IFR, and the generation time), however, the method allows back-calculating infection rates from observed death rates. Hence, in addition to the dynamic of mortality, we also reconstructed the dynamic of infection and, as a validation, compared it to external epidemiological data-like cross-sectional seroprevalence estimates-when available.

The following five parameters were estimated from the data: 1. The basic reproduction number, R 0 . According to a previous meta-analysis (Alimohamadi, Taghdir & Sepandi, 2020) , this parameter was searched in the interval 1-10. 2. The impact of non-pharmaceutical control measures, b. The lower bound of the search interval of this parameter was fixed to 0.5, such that the maximal value of the stringency index (s = 100) corresponded to a minimal reduction of SARS-CoV-2 transmission of 50% (Flaxman et al., 2020) .

3. The impact of influenza on SARS-CoV-2 transmission, β F (search interval: R). 4. The initial number of individuals exposed to SARS-CoV-2, E 1 (0) (search interval: 0-10 4 ). 5. The over-dispersion in death reporting, k D (search interval: R + ).

A summary list of fixed and estimated model parameters is presented in Table 1 . All parameters were transformed to be estimated on the real line, using a log transformation for positive parameters and the extended logistic function f (θ ) = log θ −a b−θ for parameters constrained in the interval [a,b] . The maximum iterated filtering algorithm (MIF2, Ionides et al., 2015) , implemented in the R (version 3.6.3) package pomp (King, Nguyen & Ionides, 2016; R Core Team, 2020 ) (version 2.7), was used to estimate model parameters. The R checkpoint package was used to freeze all the packages' version at the date of April 3, 2020 (De Vries, 2020) . The estimation was completed in several steps, starting with trajectory matching to identify good starting parameters for MIF2, followed by 100 independent runs of MIF2 to locate the maximum likelihood estimate (MLE). Each MIF run had 150 iterations with 5,000 particles, geometric cooling, and a random walk standard deviation of 0.1 for the initial condition E 1 (0) and of 0.02 for the other parameters. The log-likelihood of every parameter set was calculated as the log of the mean likelihood of 5 replicate particle filters, each with 20,000 particles. The profile likelihood was calculated to verify the convergence of MIF2 and to derive an approximate 95% confidence interval for the parameter β F (Raue et al., 2009 )-the parameter of key interest in our study. For the other parameters, a parametric bootstrap was used to calculate approximate 95% confidence intervals, by re-estimating the parameters for each of 200 synthetic datasets simulated at the MLE (Domenech de Celles et al., 2018; Domenech de Celles et al., 2019) . Compared with the profile likelihood, the parametric bootstrap requires less computation and was found to perform well in previous applications (Domenech deCelles et al., 2018; Domenech de Celles et al., 2019) .

To verify the robustness of the parameter estimates, we conducted six sensitivity analyses (further detailed in the Supplementary Results). First, we estimated an extended model in which the reduction of SARS-CoV-2 transmission was allowed to scale non-linearly with the stringency index. Second, to test the possible presence of other variables confounded with influenza, we estimated a model that included an exponential trend in the transmission rate of SARS-CoV-2. Finally, we varied the fixed value of 3 parameters and re-estimated the parameters of the base model. Specifically, we tested two alternative values of the average generation time (E(T g ) = 5 days and E(T g ) = 7.5 days), one alternative value of the infection fatality ratio (µ = 0.005), and one alternative value of the average onset-to-death period ( 1 κ = 13 days).

As shown in Fig. 1A , the number and the intensity of control measures against COVID-19 gradually increased from January until the nationwide lockdown in March, before a k D Over-dispersion in death reporting R + E 1 (0) Initial number exposed to SARS-CoV-2 0-10 4 Initial condition relaxation from May 2020. During this time period, the epidemic of influenza started in January 2020 and ended in March 2020, with a peak during February in each country ( Fig. 1B and Fig. S1 ). Despite correlations between some parameters (in particular the reproduction number and the impact of control measures, see Fig. S4 ), all parameters were identifiable in each country (Table 2 ). Parameter estimates indicated that, during the period of co-circulation, influenza was associated with an average 1.8-3.4-fold (uncertainty range across countries: 1.1 to 5.0) population-level increase in SARS-CoV-2 transmission ( Table 2 and Fig. S3 ). After controlling for the impact of influenza, our estimates of the basic reproduction number (R 0 ) ranged from 1.2 (in Italy) to 3.4 (in Belgium). Although the increased transmission associated with influenza early during the SARS-CoV-2 epidemic explained the data significantly better (Table 2) , a model without influenza led to higher R 0 estimates (range 2.4-5.2, Fig. 2A) , consistent with those of a previous study (Flaxman et al., 2020) . Also in line with Flaxman et al. (2020) , we found consistent evidence for a marked impact of non-pharmaceutical control measures (Table 2) , which were associated with a decrease in SARS-CoV-2 transmission below the reproduction threshold from mid-March to June 2020 ( Fig. 2A) . Notes. SE, standard error, calculated using 5 replicate particle filters, each with 20,000 particles, at the maximum likelihood estimate.

Visual inspection of simulations suggested that our model correctly captured the dynamics of COVID-19 mortality in every country (Fig. 2B) . A more detailed model-data comparison of summary statistics Wood (2010) confirmed that our model accurately reproduced the peak time, the peak number and the total number of deaths, and the death growth exponent (Maier & Brockmann, 2020) , except in Italy and Spain where the latter statistic was systematically under-estimated (Fig. S5) . Our model-based estimates of the total proportion of individuals infected with SARS-CoV-2 (as of 4 May 2020, Table 2 ) were also comparable with those of a previous modeling study (Flaxman et al., 2020) and of a seroprevalence study conducted in early May in Spain (approximate seroprevalence estimate of 5% ). Hence, our model appeared to precisely recapitulate the epidemiology of SARS-CoV-2 morbidity and mortality over a period of ∼4 months.

The results of the sensitivity analyses are presented in Tables S1-S3. We found little statistical evidence that the model with non-linear scaling of the stringency index outperformed the model with simple linear scaling in any country ( logL ∈ [0.0,1.1], likelihood ratio test P-value P ∈ [0.14,1.00], Table S1 ). Of note, although our estimates of the impact of influenza varied little, the additional estimated parameter resulted in higher parametric uncertainty, particularly in Norway where the approximate 95% CI embraced the null value. Similarly, the model with an unexplained exponential trend in transmission did not substantially improve model fit ( logL ∈ [0.0,3.2], likelihood ratio test P-value P ∈ [0.01,1.00], Table S2 ). Despite higher parametric uncertainty caused by the estimation of the trend, our results regarding the impact of influenza remained robust, but once again the approximate confidence interval embraced the null value in Norway. In addition, we found that the parameter estimates varied little when testing alternative hypotheses about the fixed value of the average generation time, of the onset-to-death time, and of the infection fatality ratio (Table S3) . Finally, because our R 0 estimates were lower than previous estimates in Italy and Spain, we tested an alternative model with R 0 fixed to 2.5 in each country. This model led to broadly similar conclusions, although the estimated impact of 

The main goal of this study was to test the hypothesis that influenza impacted the epidemiological dynamics of SARS-CoV-2, building on previous experimental evidence of a positive interaction between the two viruses (Smith et al., 2020; Ziegler et al., 2020; Bai et al., 2021) . To do so, we developed a semi-mechanistic, population-based model of SARS-CoV-2 transmission and of COVID-19 mortality, which simultaneously incorporated the impact of non-pharmaceutical control measures and of influenza. Using likelihood-based statistical inference techniques, we confronted this model with mortality incidence data in four European countries to systematically test a range of assumptions about the possible impact of influenza and of control measures. In keeping with previous studies (Flaxman et al., 2020) , we found robust and consistent evidence that control measures markedly reduced the transmission of SARS-CoV-2. In addition, we also found consistent evidence suggesting that co-circulation of influenza transiently facilitated the transmission of SARS-CoV-2 early in the epidemic in Europe. Our study has a number of important limitations. First, as in other studies (Flaxman et al., 2020; Li et al., 2020; Kucharski et al., 2020) and because of a lack of appropriate age-specific data (for example on the temporal changes in the contact matrix and in the incidence of influenza), our model was not age-structured, even though many aspects of COVID-19 and of influenza epidemiology-like disease severity and lethality-vary markedly with age (Verity et al., 2020) . The susceptibility to SARS-CoV-2 infection was also found to increase with age (Davies et al., 2020) , a finding potentially explained by lower baseline expression of the ACE2 receptor in children (Bunyavanich, Do & Vicencio, 2020) . A testable prediction of our model, therefore, is that influenza should be associated with a transient increase in susceptibility to SARS-CoV-2 infection, commensurate with the variations of influenza incidence over age. Second, we modeled the impact of nonpharmaceutical control measures using a simple, linear function scaling the stringency index to the reduction of SARS-CoV-2 transmission. Even though this simple hypothesis provided a more parsimonious fit, that result may be specific to Europe, where control measures gradually increased in number and in intensity (Fig. 1A) . In general, the association is likely non-linear (e.g., if a high-impact intervention like a lockdown is implemented early on), and we therefore recommend testing a variety of scaling functions. More generally, although our model builds on a previously validated method to estimate the time-varying reproduction number (Flaxman et al., 2020) , we acknowledge that the stringency index may not fully capture temporal variations in SARS-CoV-2 transmission, in particular behavioral changes outside of what was mandated by governments. Even though we tested a model with an unexplained trend in transmission (Table S2) , more complex temporal functions may be required to fully capture such changes. Third, we did not incorporate climate into our model, even though, as for other respiratory viruses, environmental variables like temperature and humidity may affect the transmission of SARS-CoV-2. According to previous studies conducted in a variety of locations worldwide, however, the impact of weather on the SARS-CoV-2 epidemic appears to have been modest, at least during the first wave in early 2020 (Jüni et al., 2020; Gaudart et al., 2021; Sehra et al., 2020) . These findings are also consistent with epidemiological theory, which predicts that, because of lack of population immunity, the initial pandemic trajectory may be relatively insensitive to climate (Baker et al., 2020) . Fourth, we did not specifically model fully asymptomatic cases, which may represent a large fraction of SARS-CoV-2 infections (Li et al., 2020) . The omission of asymptomatic infections may lead to biased R 0 estimates if their duration significantly differs from that of symptomatic infections (Park et al., 2020) . A previous study, however, estimated that the duration of both types of infection is comparable (Li et al., 2020) , such that our estimates should be robust in more complex model structures. Finally, we assessed only the impact of influenza, because of its high prevalence and period of overlap with SARS-CoV-2 in early 2020 in Europe and of the availability of high-quality data (Kramer & Shaman, 2019) . Nevertheless, other respiratory viruses, like RSV and rhinoviruses (Dee et al., 2021) , may also interact with SARS-CoV-2 and could be considered.

Acknowledging these limitations, our model makes at least two other predictions that could be tested to provide confirmatory evidence. First, even though our results did not allow to distinguish between higher transmissibility or higher susceptibility in individuals co-infected with influenza and SARS-CoV-2, previous experimental work suggests that the latter mechanism may operate, as a result of up-regulation of the ACE2 receptor caused by influenza infection (Smith et al., 2020; Ziegler et al., 2020) . Hence, we predict that a recent influenza infection should be an independent risk factor for subsequent SARS-CoV-2 infection. Estimates of the frequency of co-detection of influenza and SARS-CoV-2 by polymerase chain reaction (PCR) testing in nasopharyngeal swabs were highly variable in previous studies (range 0-60% (Thindwa et al., 2020; Ozaras et al., 2020) ). Although the marked seasonality of influenza in temperate regions likely explains in part the low frequency found in some studies (Ozaras et al., 2020) , we propose that differences in the natural history of influenza and SARS-CoV-2 infections also lead to a systematic under-estimation of co-infection. Specifically, because the incubation period of SARS-CoV-2 infection (estimated to average 5.7 days (Khalili et al., 2020) ) exceeds that of influenza (A, 1.4 days or B, 0.6 days (Lessler et al., 2009)) , it is likely that, by the time SARS-CoV-2 infection becomes detectable, influenza no longer is. To make that statement more precise, we calculated the probability of detectability of a co-infection, with influenza first then SARS-CoV-2 (Table S4 ). Assuming that influenza is detectable by PCR up to 4-5 days after (Carrat et al., 2008) , and SARS-CoV-2 from 2-4 days before (Tindale et al., 2020) , symptom onset, we find that a large fraction (30-50%) of co-infections may not be detectable at all. These results may help explain the low frequency of co-detection found in some studies (Kim et al., 2020) , and suggest that the time window of co-detectability may be too short to adequately infer the association between influenza and SARS-CoV-2 using PCR testing. Serological studies comparing the prevalence of antibodies against influenza in SARS-CoV-2 cases and non-cases may therefore be required to test the prediction that influenza is a risk factor for SARS-CoV-2 infection. Second, we predict that, all else being equal, individuals vaccinated against influenza should be at lower risk of SARS-CoV-2 infection than those unvaccinated. The findings of a negative association between influenza vaccine coverage and COVID-19 mortality in ecological studies (in Italy (Marín-Hernández, Schwartz & Nixon, 2020) and in other countries (Arokiaraj, 2020)) and of a lower risk of SARS-CoV-2 infection in influenza vaccinees in some individual-level epidemiological studies (reviewed in Riccio et al., 2020) are consistent with our prediction, but further epidemiological investigations are needed. Importantly, our results can explain these findings as the direct effect of influenza vaccines on influenza infection, instead of indirect effects on non-influenza pathogens (e.g., as a result of trained immunity) (Salem & El-Hennawy, 2020) .

With the likely prospect of COVID-19 becoming endemic, the potential interactions of SARS-CoV-2 with other respiratory pathogens-in particular respiratory viruses-may become a key public health issue. In this context, our results suggest that influenza could facilitate the circulation of SARS-CoV-2 and therefore increase the burden of COVID-19. As outlined above, these results are consistent with several lines of experimental (Bai et al., 2021; Smith et al., 2020; Ziegler et al., 2020) and epidemiological (Riccio et al., 2020; Arokiaraj, 2020; Marín-Hernández, Schwartz & Nixon, 2020) evidence. We note, however, that a previous study proposed that influenza and SARS-CoV-2 have competitive interactions (Pinky & Dobrovolny, 2020) . Specifically, using a within-host model of viral replication Pinky and Dobrovolny found that the low growth rate of SARS-CoV-2 may result in limited access to target cells and therefore suppression by other respiratory viruses (Pinky & Dobrovolny, 2020) . Although this mechanism may be generally relevant for respiratory viruses, the ability of influenza to up-regulate ACE2 (Smith et al., 2020; Ziegler et al., 2020 )-a feature not included in the within-host model (Pinky & Dobrovolny, 2020 )-could counteract this mechanism and explain the increased infectivity of SARS-CoV-2 found experimentally (Bai et al., 2021) . Of note, in keeping with Pinky & Dobrovolny (2020) , another experimental study found evidence that the interferon response caused competitive interactions between rhinoviruses and SARS-CoV-2 (Dee et al., 2021) . Hence, these different biological mechanisms suggest that every respiratory virus may interact with SARS-CoV-2 in a highly specific way, with influenza being unique in its ability to up-regulate ACE2 and to increase SARS-CoV-2 infectivity (Bai et al., 2021) .

In conclusion, our results suggest that influenza virus infection could have increased the transmission of SARS-CoV-2 and facilitated its spread during the early 2020 epidemic of COVID-19 in Europe. Hence, an increase in the uptake of influenza vaccines may be called for, not only to reduce hospitalizations due to influenza infections (Ozaras et al., 2020; Paget et al., 2020) , but also to reduce their downstream impact on SARS-CoV-2 transmission and on COVID-19 mortality. More generally, taking into account the microbial environment

Estimate of the basic reproduction number for COVID-19: a systematic review and meta-analysis

CMMID COVID-19 working group, Eggo RM. 2020. Age-dependent effects in the transmission and control of COVID-19 epidemics

Human rhinovirus infection blocks SARS-CoV-2 replication within the respiratory epithelium: implications for COVID-19 epidemiology

The association between influenza vaccination and the risk of SARS-CoV-2 infection, severe illness, and death: a systematic review of the literature

Microsoft, Checkpoint: install packages from snapshots on the checkpoint server for reproducibility

COVID-19 italia-monitoraggio della situazione

Unraveling the seasonal epidemiology of pneumococcus

The impact of past vaccination coverage and immunity on pertussis resurgence

Estimating the overdispersion in COVID-19 transmission using outbreak sizes outside China

Geographic distribution of COVID-19 cases worldwide

Estimating the effects of non-pharmaceutical interventions on COVID-19 in Europe

Factors associated with the spatial heterogeneity of the first wave of COVID-19 in France: a nationwide geo-epidemiological study

Centre for Mathematical Modelling of Infectious Diseases COVID-19 working group. 2020. Early dynamics of transmission and control of COVID-19: a mathematical modelling study

Incubation periods of acute respiratory viral infections: a systematic review

Substantial undocumented infection facilitates the rapid dissemination of novel coronavirus (SARS-CoV2)

Destabilization of epidemic models with the inclusion of realistic distributions of infectious periods

Superspreading and the effect of individual variation on disease emergence

Effective containment explains subexponential growth in recent confirmed COVID-19 cases in China

The impact of pandemic influenza A (H1N1) 2009 on the circulation of respiratory viruses

Epidemiological evidence for association between higher influenza vaccine uptake in the elderly and lower COVID-19 deaths in Italy

Coinfection can trigger multiple pandemic waves

Predominant role of bacterial pneumonia as a cause of death in pandemic influenza: implications for pandemic influenza preparedness

Age-specific mortality and immunity patterns of SARS-CoV-2

Influenza interaction with cocirculating pathogens and its impact on surveillance, pathogenesis, and epidemic profile: a key role for mathematical modelling

An open call for influenza vaccination pending the new wave of COVID-19

The impact of influenza vaccination on the COVID-19 pandemic? Evidence and lessons for public health policies

Infection fatality risk for SARS-CoV-2 in community dwelling population of Spain: nationwide seroepidemiological study

SARS-CoV-2 coinfections: could influenza and the common cold be beneficial?

Prevalence of SARS-CoV-2 in Spain (ENE-COVID): a nationwide, population-based seroepidemiological study

R: a language and environment for statistical computing

Vienna: R Foundation for Statistical Computing

Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood

The possible beneficial adjuvant effect of influenza vaccine to minimize the severity of COVID-19

Estimating the burden of SARS-CoV-2 in France

Maximum daily temperature, precipitation, ultraviolet light, and rates of transmission of severe acute respiratory syndrome coronavirus 2 in the United States

Cigarette smoke exposure and inflammatory signaling increase the expression of the SARS-CoV-2 receptor ACE2 in the respiratory tract

A note on generation times in epidemic models

Use of seasonal influenza and pneumococcal polysaccharide vaccines in older adults to reduce COVID-19 mortality

Evidence for transmission of COVID-19 prior to symptom onset

Estimates of the severity of coronavirus disease 2019: a model-based analysis

Appropriate models for the management of infectious diseases

Factors associated with COVID-19-related death using OpenSAFELY

Statistical inference for noisy nonlinear ecological dynamic systems

China Novel Coronavirus Investigating and Research Team. 2020. A Novel Coronavirus from Patients with Pneumonia in China

SARS-CoV-2 receptor ACE2 is an interferon-stimulated gene in human airway epithelial cells and is detected in specific cell subsets across tissues

of SARS-CoV-2 may be essential, not only to better understand its epidemiology, but also to enhance current and future infection control strategies.

No specific funding was used for this study.

Matthieu Domenech de Cellès received postdoctoral funding (2017-2019) from Pfizer and consulting fees from GSK. Jean-Sebastien Casalegno declares no competing interests. Bruno Lina is chair of the ISC for the Global Influenza Surveillance Network, and co-chair of the Global Influenza and RSV initiative; he is also a member of the French COVID-19 Scientific Committee (no personal income for all these activities). Lulla Opatowski has received consulting fees from WHO for work on antimicrobial resistance and funding from Pfizer (2017-2019).

• Matthieu Domenech de Cellès conceived and designed the experiments, performed the experiments, analyzed the data, prepared figures and/or tables, authored or reviewed drafts of the paper, and approved the final draft.• Jean-Sebastien Casalegno, Bruno Lina and Lulla Opatowski conceived and designed the experiments, authored or reviewed drafts of the paper, and approved the final draft.

The following information was supplied regarding data availability:The data used for this analysis are freely available from the databases cited in the references and have been compiled into a database available from Edmond, the Open Data Repository of the Max Planck Society: https://dx.doi.org/10.17617/3.7r. All R programming codes are also available from this repository.

Supplemental information for this article can be found online at http://dx.doi.org/10.7717/ peerj.12566#supplemental-information.