## Could short-lasting non-specific immunity explain seasonal spacing of epidemic diseases?

Gorm Gruner Jensen, 1, ∗ Florian Uekermann, 2, ∗ Kim Sneppen, 2 and Lone Simonsen 3

1 Institute for Theoretical Physics, University of Bremen, D-28359 Bremen, Germany 2 Niels Bohr Institute, University of Copenhagen, Copenhagen 2100-DK, Denmark 3 Department of Science and Environment, Roskilde University, DK-4000 Roskilde, Denmark

Common respiratory viruses cause seasonal epidemics in sequential patterns. The underlying mechanisms for this pattern have been debated for some time. For influenza, contenders include temperature, humidity, and vitamin D levels. While such seasonal drivers may be sufficient to explain midwinter peaking it is unclear if other respiratory viruses peaking in spring, summer, or fall have individual environmental drivers. Here we present a dynamic model of interacting diseases, in an effort to explain observed seasonal patterns without requiring multiple seasonal drivers. Our model extends the classical SIRS-model to include multiple diseases that interact via a short-lasting, nonspecific immunity component. This temporal protection is triggered by exposure to any epidemic disease and lasts only a couple of weeks. We show that the inhibiting disease-interaction allows recurrent epidemic behaviour without any seasonal driving. In the presence of a single seasonal driver our model parsimoniously predicts epidemic patterns such as sequentially occurring epidemic diseases. We also present a two-disease simulation reproducing multiple features observed in time series of parainfluenza strains PIV-3 and biennial PIV-1 epidemic patterns. Complex epidemic patterns of seasonal diseases may be explained by non-specific innate or cloub cell (T-cell) mediated immune responses that result in interaction between unrelated respiratory epidemic diseases.

## KEYWORDS

Seasonal epidemics; Seasonal ordering; Non-specific immunity; Disease interaction; SIRS model;

## I. INTRODUCTION

Mysteriously, each year respiratory viruses cause epidemics that peak in predictable months during spring, summer, fall and winter [1-6]. The fact that the Northern and Southern hemisphere displays mirror image patterns of timing suggest that epidemic timing is determined by seasonal environmental driver(s). This has been known probably since the beginning of virology in the 1930s, but the phenomenon is still in search of an explanation. Various authors have suggested that winter-seasonality in for example influenza is explained by a seasonal driver such as cold weather, low humidity or low vitamin D levels [1, 2, 7-10]. However, in order for such mechanisms to explain spring, summer, and autumn seasonality one would need to imagine a different seasonal driver for each season - perhaps for each virus - not an appealing or parsimonious hypothesis. Rather, the repeated annual patterns suggest that somehow epidemic diseases can repel each other and await their turn. Such a mechanism of epidemic interference has been proposed by Rohani [11, 14, 16] who demonstrated that by introducing short-lasting, non-specific immunity into classical dynamic childhood disease models one can find parameters for which one epidemic would exclude another. Yet most mathematical models of disease-

∗ G. G. J. and F. U. contributed equally to this work.

disease interactions published to this date explore the effects of Ab-mediated specific, life long cross immunity [12, 13, 15, 17]. When short-lasting immune mechanisms are considered it is typically in a single disease context to explain phenomena like multi-wave diseases [18-22] or to provide an evolutionary bottleneck for influenza evolution [23]. Thus, models have yet to explored the possibilities of inter-epidemic dynamics in multiple-pathogen models that include terms of short-lasting non-specific immunity. Here we address this intriguing question anew, proposing a dynamic disease model form that can theoretically address this type of inter-species interaction. Our model work is motivated by the ever increasing evidence for short-lasting acquired non-specific immunity [22, 24-27].Prominently, the mechanism has been directly demonstrated in controlled laboratory experiments with influenza A infections that prevented subsequent infection with influenza B for 6 to 12 weeks in mice [27].

We explore theoretically whether it is possible to explain the phenomenon of recurrent multiple epidemics with predictable (and different) seasonality by introducing short-term non-specific immunity into the classical type of SIRS modelling and choosing parameter values consistent with known mechanisms such as innate immunity. This way we allow for individuals infected with a disease to be temporarily unavailable for infection with a second disease, and then re-enter the susceptible pool for the second disease after a few weeks. To test the performance of our model against observed epidemic patterns, we apply the model form to real-life epidemiological time series data of PIV-1 and PIV-3 in the US. We find that our model can reproduce the observed signature patterns of seasonal diseases assuming only a single environmental driver and one type of short-lasting non-specific immunity.

Before introducing the model in detail we emphasise

R2

ase 2

Figure 1. Model from the perspective of individual hosts. Each host is in either the I , R , or S state with respect to each disease, and either has non-specific immunity ( G on ) or not ( G off ). The green arrows represent exposure of host 1 to disease 1 due to contact with the infected host 2. Infections ( S → I ) require exposure to the disease, but are prevented if the the host has non-specific immunity G on (as indicated by the red links).

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00925v1-figure-aac8f234-92e3-40e2-93a6-aedcc2d06618.jpg)

that there are multiple biological motivations for short non-specific immunity. While it may represent a universally active innate immunity [24], any other less specific short-term immune defences that affect a class of diseases [27, 28] are valid motivations as well. We will not explicitly distinguish between these motivations in the following, as our model is general enough to represent any of them. As a result, it must be decided on an individual basis whether sub-types or strains of viruses should be considered distinct diseases in the context of this model or not. The crucial criteria are that diseases should not confer significant long-term cross-immunity to each other and allow for reinfection. For example, in the case of influenza A, we would consider the different sub-types (H3N2,H1N1) to be distinct diseases, since infection with one does not typically confer long-term immunity against another. The strains of a sub-type however would be considered the same disease here, to allow reinfection every couple of years. Under that representation influenza 'waning' specific immunity become an abstraction of antigenic drift and cluster transitions. [29] With respect to PIV (parainfluenza virus) and RSV (respiratory syncytial virus) on the other hand, the waning immunity represents an actual loss of effective specific immunity in the host. [30, 31]

## II. MODEL AND METHODS

Our model is based on the traditional SIRS -model[32] with a homogeneous, fully mixed host population. We extend this model to multiple diseases and add a new feature, the effect of a non-specific immune protection, which is short-lived. We think of this additional kind of immunity as a way of modelling either an increased efficiency of the innate immune system observed by [27] or as short-lived cross-immunisation between similar, but not identical, diseases. All events, like contact between hosts or recovery from a disease, are modelled on the basis of transition rates. In agent-based simulations with hostpopulations of finite size this leads to stochastic dynamics, while the behaviour for very large host-populations can be approximated using differential equations. We will describe the model from the perspective of individual hosts, which is shown schematically in figure 1.

SIRS: The traditional SIRS -model can be summarised as follows: Each host can be S usceptible, I nfectious, or R emoved. The transmission of diseases is modelled using a contact rate. If the contacted host is susceptible, it becomes infectious ( S → I ). During one period of being infectious, a host will on average contact R 0 randomly chosen hosts. This is commonly referred to as the basic reproduction number. The average time a hosts stays infectious until entering the removed state ( I → R ) is T I . Similarly, a removed host becomes susceptible again ( R → S ) after an average time of T R .

Multiple diseases: To describe multiple diseases spreading in the same host population, we track whether the host is S usceptible, I nfectious, or R emoved with respect to each disease independently (see disease 1 and 2 in Fig. 1). The R state therefore represents a disease specific immunity.

We then add a non-specific immunity, which is either active ( G on ), or inactive ( G off ). If it is active, the host is protected against infection with any disease, even if it is in the S state with respect to some of them. The nonspecific immunity is activated ( G off → G on ) every time a host comes into contact with any of the diseases, also if infection is prevented by specific immunity (host is in R state with respect to the disease). This assumption allows maximal interaction between diseases. The average time the non-specific immunity stays active is T G ( G for G eneral Immunity).

Since the assumption of short-term but 100% effective non-specific immunity is not necessarily realistic, we conducted simulations with non-specific immunity that is twice as long, but not fully protective (50%). We achieved both qualitatively and quantitatively similar results (data not shown).

Model parameters: The model has 4 parameters, three of which are independent. Throughout the paper we keep as many parameters as possible fixed to the same value.

R 0 : Basic reproduction number. Mean number of transmission attempts over the course of one infection. Expected number of secondary cases in a completely susceptible population. R 0 = 2 5 in case . studies without seasonal variation. Case studies with a seasonal driver use R 0 ( ) = 2+ t c · sin(2 πt/τ ), with c = 0 1 or . c = 0 2 and . τ = 1year.

T I : Mean duration of infectiousness. Since scaling all durations T X by a common factor has no dynamic effect, we set T I = 1 without loss of generality. For common respiratory diseases, T I = 1 represents ≈ 6

days. All simulations are conducted with T I = 1, without exception. All other durations are specified in multiples of T I .

T R : Mean duration of disease specific immunity. Unless stated otherwise T R = 100.

T G : Mean duration of short non-specific immunity T G ≤ 6. See individual figure for specific values.

R 0 , T I and T R can be disease specific parameters. Unless stated otherwise, they are chosen to be identical for all diseases.

The mean duration of specific immunity is rather short in most simulations (few years). This reflects the fact that immunity is often incomplete in viruses such as PIV[30] and RSV[31], which allows for frequent reinfection. In influenza on the other hand, the effective duration of specific immunity is shortened by the antigenic evolution of influenza[29].

Simulation details: All simulations discussed in this paper were initialised with an almost exclusively susceptible population. Stochastic agent-based simulations with a finite number of hosts were performed using a Gillespie algorithm. The population consists of 500 000 susceptible hosts and a single 'exterior host', who is infectious with all diseases and never recovers. The 'exterior host' represents external sources of infections, which has the side-effect of preventing stochastic extinction.

Similarly, the deterministic simulations contain an 'external reservoir' of always infectious hosts which is 500 000 times smaller than the simulated population. The external reservoir dramatically reduces transient time, i.e. the time before a stable pattern emerges. In simulations with a seasonal driver the external reservoir is 200 000 times smaller than the simulated population. This increase of external infection severely reduces the transient time before the epidemic limit-cycle is reached.

The deterministic rate equations are solved using the Euler method with a time-step of 0 001 . T I . To break the symmetry between diseases with identical dynamical parameters, we initialise the deterministic simulations with a one thousandth of the population being infectious with disease one, and two thousandths of the population being infectious with disease two, etc.

We have not systematically explored the robustness of the dynamics to initial conditions. However, we also have not observed sensitivity to the choice of initial-conditions, except for a small parameter region of bi-stability marked in Fig. 2c.

## III. RESULTS

First, we investigate the simplest case where interactions can be observed: Two diseases with identical dynamical parameters, without changing environmental factors. Fig. 2a&amp;b shows two case studies with different durations of average of non-specific immunity ( T G ), for which we observe qualitatively different behaviour. In the case of very short non-specific immunity ( T G = 0 1), both . diseases approach an endemic equilibrium state as shown in Fig. 2a. This resembles the well known behaviour without non-specific immunity (traditional SIRS ). The other case, shown in Fig. 2b, considers a longer duration of non-specific immunity ( T G = 4, corresponding to about 24 days for a disease with a 6-day infectious period). Here we observe a repeated pattern of distinct epidemics, alternating between the two diseases. The agentbased and deterministic (infinite population size) simulations exhibit the same qualitative behaviour. These examples show that introduction of non-specific immunity can lead to recurrent alternating epidemics, a feature that unmodified deterministic SIRS models do not exhibit[33]. These case studies show that the non-specific immunity needs to exceed a certain average duration to maintain recurrent epidemics.

In Fig. 2c we further investigate the phenomenon by measuring the amplitude of the recurring epidemics for specific choices of R 0 and duration of non-specific immunity T G . Since realistic durations for the non-specific immunity are short (2 -4 weeks), we only show the behaviour for T G &lt; 15, which corresponds to up to 3 months, if the infectious period is 6 days long. In this range, recurrent epidemics are observed for most reasonable values of R 0 and T G . More specifically, recurrent epidemics are only absent if the interaction between the two diseases is very weak (short non-specific immunity time) or R 0 is close to 1.

The time series in Fig. 2a&amp;b show the transient behaviour of simulations. Much longer simulations, such as those underlying Fig. 2c show that the behaviour indeed converges to stable periodic attractors.

## III.1. Recurring alternating epidemics

For a given parameter set, we can numerically find the endemic fixed point and then perform a linear stability analysis. We have done this on a fine grid over the parameter spaces shown in Fig. 2c &amp; 3c and found supercritical Hopf-bifurcations, where the stable (endemic) fixed-point becomes unstable and emits a stable (epidemic) limitcycle. The bifurcation lines are shown as dashed white lines in Fig. 2C &amp; 3C. The amplitudes measured in deterministic simulations match well with the calculated bifurcation lines.

## III.2. Diseases with different dynamical parameters

In the previous section we established the emergence of recurrent epidemics for two diseases for a wide range of dynamical parameters, which were set to be identical for both diseases. A natural concern is, whether this behaviour is robust to variations of dynamical parameters between the two diseases. Next we explore the dynamics of two diseases with differing parameters for transmissibility ( R 0 changes, while T I = 1 remains unchanged)

Infectious

Figure 2. Dynamics with two diseases with mean specific immunity duration of T R = 100 T I . (a&amp;b) : Case studies with R 0 = 2 5. . Continuous lines show stochastic agent-based simulations with 500000 hosts. Dotted lines show one of the diseases in deterministic simulations. (a) : Short non-specific immunity T G = 0 1 . T I , leads to endemic equilibrium state. (b) : Longer non-specific immunity T G = 4 T I , leads to recurrent alternating epidemics. (c) : Spikiness of epidemic behaviour. For each data point the deterministic simulation of two diseases is run for a period of 11000 T I . The maximum and minimum of infectious hosts ( I max , I min ) during the last 1000 T I is measured and log 10 I max I min is plotted. The white dashed line was obtained by linear stability-analysis and marks the bifurcation line between stable endemic fixed point and epidemic limit-cycles. The black dotted line encloses an area with two stable epidemic limit-cycles (the plotted color corresponds to the one with larger amplitude).

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00925v1-figure-2f2f6e12-5d18-437c-973a-8c912e3d7ab5.jpg)

and duration of specific immunity ( T R ). We approach this question by varying the disease specific R 0 ,d and the duration of specific immunity T R,d of the two disease in opposite directions:

$$R _ { 0, d = 1 } & = R _ { 0 } + \Delta R _ { 0 } \quad \quad R _ { 0, d = 2 } = R _ { 0 } - \Delta R _ { 0 } \left ( 1 \right ) \quad \text{the s\nmid} \\ T _ { R, d = 1 } & = T _ { R } + \Delta T _ { R } \quad \quad T _ { R, d = 2 } = T _ { R } - \Delta T _ { R } \left ( 2 \right ) \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \ \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \text{susta} \quad \end{$$

We choose R 0 = 2 5 and . T R = 100 T I as base values and set the average duration of non-specific immunity to T G = 6 T I . As in the previous section, we will first examine two case studies, which are shown in Fig. 3a&amp;b. The exact choice of parameters is shown in Fig. 3c. In the first case study R 0 and T R are both increased for one disease and decreased for the other, which leads to a dynamics of alternating recurrent epidemics, similar to the pattern observed for diseases with identical parameters. In the second case study the two parameters are varied in opposite directions, R 0 is increased, while T R is decreased for one disease and vice-versa for the other. The epidemic pattern still shows recurring epidemics, however, the two disease exhibit a different number of peaks per period and less pronounced peaks.

Infectious

Figure 3. Dynamics with two diseases with different dynamical parameters. The transmissibility R 0 and specific immunity duration T R differ between diseases by ± ∆ R 0 and ± ∆ T 0 from their respective mean R 0 = 2 5 and . T R = 100 T I . T G = 6 T I remain unchanged, and T I is the same for both diseases. (a&amp;b) : Case studies. Stochastic simulations are shown in solid lines. Dashed line shows the deterministic results for one of the diseases in b. ∆ T 0 and ∆ R 0 are marked in c in the corresponding line color. (c) : Spikiness of epidemics. See Fig. 2 for details on measurement.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00925v1-figure-599c4185-7479-44c1-8f30-016b47f275a3.jpg)

We then explore this phenomenon in more detail (figure 3c). In general larger R 0 increases the frequency of recurrent epidemics in our model, while larger T R reduce the frequency. The diagonal region with high amplitude (red) in the centre of Fig. 3c signifies parameter variations, with little change in peak frequency, where the differences in R 0 and T R compensate each other. This produces alternating recurrent epidemics as seen in 3a. Parallel to this central region are two secondary areas with high amplitudes, which represent parameter sets, where one peaks twice as often as the other, as seen in the second case study. If the periods don't match, the sustained alternating epidemics are replaced by a more varied temporal pattern with lower amplitudes.

We find that the epidemic behaviour is robust to moderate differences in both R 0 and T R , especially when the variations are such that the natural periods of recurrence of the two diseases are similar or approximately multiples of each other.

## III.3. Multiple Diseases

After exploring the case of two diseases in detail, we will briefly show examples with more than two diseases. The model is defined for an arbitrary number of diseases, so it's straight forward to extend our investigations. The number of possible host-states however, grows exponen-

Infectious

Figure 4. More than 2 diseases in deterministic simulations with T G = 6 T I . The basic reproduction number is R 0 = 2 5 except for panel (d), where . R 0 = 5. The mean specific immunity time is T R = 100 T I .

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00925v1-figure-5187ae3d-5d75-444a-af44-1dcd9ee32c00.jpg)

tially with the number of diseases, which quickly makes the simulations unfeasible in the deterministic limit of infinitely many agents. 'Individual host'-based simulations could be a way around this issue, but it would not be possible to study the model without noise. Deterministic simulations for up to 5 diseases are shown in Fig 4. With increasing number of diseases the temporal pattern becomes more varied. In general only systems with a limited number of diseases can sustain alternating recurrent epidemics.

We observe that larger values of both R 0 and T R allow for narrower epidemic spikes, which increases the number of supported diseases without breaking the pattern of alternating epidemics. Narrower epidemic spikes are necessary to sustain many alternating recurrent epidemics, because the non-specific immunity suppresses simultaneous outbreaks of different diseases, thus limiting the fraction of time available to the individual disease.

## III.4. Seasonality

A common method for incorporating seasonal factors in SIRS -models is by making R 0 time-dependent. Consequently, this section will discuss an extension of our model with an annually varying basic reproduction number R 0 ( ) t = 2 + 0 1 . · sin(2 πt/τ ), where τ specifies the length of a year in multiples of T I .

Fig. 5 shows different 3-disease patterns produced by our model for different choices of the duration of nonspecific immunity T G . In this context we can think of T G as a control-parameter for the strength of interaction between the diseases. It is worth noting that a similar set of patterns can be found by varying T R , or the meanvalue or amplitude of R 0 ( ). t We use the same dynamical parameters for all three diseases and subject them to the same seasonal driving via R 0 ( ). t Similar to the multitude of observed epidemic patterns in nature, our model produces a variety of dynamics, based on the exact parameter choices. In the following we will relate some patterns observed in nature to the different panels in Fig. 5:

In the case of concurrently spreading influenza A sub-

Infectious

Figure 5. Deterministic simulations with three diseasesystems with different strength of non-specific immunity T G . All dynamical parameters except for T G are equal in all panels and for all diseases: R 0 ( ) = 2+0 1 sin(2 t . · πt/ 120), T R = 250 T I . (a) With weak interaction (low values of T G ), the seasonal driver dominates the dynamics. (b&amp;c) With intermediate interaction-strength the seasonal driver controls the periodicity, but the epidemic peaks are separated in time. (d) Strong interaction prevents completion of distinct epidemics during every year. In combination with the seasonal driver this leads to a loss of the annual periodicity.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00925v1-figure-ef57151d-3c89-443c-a2a7-052337885b11.jpg)

types (H3N2 and H1N1) whose concurrently dominant strains have different antigenic profiles, we see that the dynamics is dominated by the seasonal driver, with few differences between disease variants. Unsurprisingly, this is reproduced in panel a, where the non-specific immunity is short, which implies weak interaction.

Some types of PIV on the other hand, tend to peak at different times than others. PIV-3 typically peaks in early summer, while PIV-1, 2 and 4 are more common in the fall. Panel B shows a similar pattern between diseases with medium duration of non-specific immunity. One disease peaks separately from the other two, which remain synchronised.

Panel C shows the behaviour with non-specific immunity that lasts long enough to spread the peaks of all three diseases over more than half a year, even though they have the same dynamical parameters and depend on the same seasonal driver. Finally, in panel D, the interaction between diseases is so strong that the epidemic behaviour does not converge to a pattern with a one year periodicity.

## III.5. Modelling observed PIV-1&amp;3 seasonal pattern

For PIV-1 and PIV-3 a particularly interesting pattern is typically observed. The activity of PIV-3 between epidemic peaks is consistently lower every second year, which coincides with the biennial outbreaks of PIV1 [4, 5]. It has been suggested that this correlation is caused by an interaction between the two diseases[4].

We found that our model can reproduce multiple features of the described dynamics, with a seasonally varying R 0 and a longer specific immunity-time for PIV-1 than for PIV-3 (see Fig. 6 for detailed information on

Figure 6. Upper panel: Seasonal variant of the deterministic model with two diseases, resembling PIV-1 and PIV-3. The basic reproduction number is varied seasonally using a sine function R 0 ( ) = 2 + 0 2 t . · sin(2 πt/ (60 T I )). The average immunity time is T G = 4 T I and the specific immunity times are T R,PIV 1 = 170 T I and T R,PIV 1 = 120 T I . Rescaled to a year, this implies an infectious period of ≈ 6 days, on average ≈ 24 days of non-specific immunity, and specific immunity that lasts on average 2 years (PIV-3) and ≈ 2 8 years (PIV-1). . The external reservoir of constantly infectious hosts is 1 / 200 000 of the population size. The annual maxima of R 0 are aligned with March 1st in the lower panel. Lower panel: Percentage of positive tests for PIV-1 and PIV-3 serotypes in the US reported to NREVSS from July 1990 to June 2004. This data was extracted from Fry et al. [4]. Due to the low image quality of the original, the estimated extraction uncertainty is ( ± 10days , ± 0 2%). .

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00925v1-figure-250fbf8b-340f-4666-86d3-68e6f020c74a.jpg)

parameters). The direct comparison between simulated and observed prevalence in Fig. 6 shows that the model reproduces many features of the observed dynamics, such as: (1) PIV-3 peaks annually, while PIV-1 peaks every second year. (2) In years with a PIV-1 epidemics, the valley of PIV-3 activity is significantly deeper. (3) The next PIV-3 epidemics is delayed by ≈ 3 months in years with a PIV-1 epidemic. (4) PIV-1 epidemics ascend for many months and end abruptly, while the individual PIV-3 peaks are more symmetrical. If the simulation is executed with zero - or very low levels of - non-specific immunity, but otherwise identical parameters, the system converges to a limit-cycle in which both diseases have epidemic outbreaks every year.

between diseases is too weak to overcome the effect of the environmental driver. This seems to indicate, that the observed effect would not occur just because of behavioural changes caused by the sickness, such as staying home from work.

## IV. DISCUSSION

Many respiratory viral diseases cause seasonal epidemics in predictable months. It is well known that epidemics of different diseases spike at different times of the year, but it is not clear why. Typical discussions of epidemic patterns point to a variety of seasonal drivers[1, 2, 7-10]. It is not clear, however, why these would have a different influence on otherwise similar respiratory pathogens causing these epidemics. Here we have demonstrated that a simple SIRS-model with a single seasonal driver can predict different seasonal patterns even for diseases with identical dynamical parameters (figure 5), provided that these diseases interact through a short-lasting non-specific immune response in their hosts.

At an organism level, innate immunity, or crossimmunity between specific viral strains, may have many different causes[22, 24-27]. Our non-specific immunity mechanism could be interpreted as any of these or even a combination. It is, however, worth noting that the effect should be short-lasting, i.e. much shorter than the effect of antibody mediated (pathogen specific) immunity, but also not too short. In figure 5a it is observed that if the average duration of the non-specific immunity is as short as the average infectious period, the interaction

This paper is not the first to make a theoretical investigation of how non-specific immunity may change the epidemic patterns of diseases. Some works have explored a similar concept using a coalescence, following the infection of one disease, in which hosts are neither infectious nor susceptible to any diseases[11, 14, 16]. However, these models are mostly inspired by phenomena related to childhood diseases, which leads to some fundamental differences in the underlying single-disease dynamics. In particular, the pathogen-specific immunity is assumed to last forever, making it necessary to explicitly keep track of birth and death processes. In addition, these models have, to the best of our knowledge, not yet been used to directly address the question of seasonal ordering.

We have found that our model is able to predict recurrent epidemic oscillations, even in the absence of seasonal driving. This is an interesting result in itself, as it is well known that epidemic limit cycles don't exist for any parameters in the classical SIRS model. In the figures 2 and 3, we have investigated these recurrent epidemic patterns using a combination of local linear stability analysis and direct simulations. The reason why non-specific immunity can break the stability of the endemic fixed point in our model is that one disease postpones the spread of others, thereby increasing the number of hosts without specific immunity to these temporarily suppressed diseases. In turn, this causes larger epidemics once the non-specific immunity wears off, thus leading to epidemic peaks, as opposed to stable levels infectious hosts. While these results are interesting, they are not particularly surprising since similar phenomena have been observed in related models with short term non-specific immunity[11, 14, 16], as well as epidemic models with life long cross immunities[12, 13, 15, 17].

Finally we have shown an example of how our model is able to reproduce a number of qualitative features in

the complex epidemic pattern of PIV-1 and PIV-3 observed in the period between 1990 and 2004 in the United States. This strengthens the hypothesis proposed by Fry et. al[4], that this pattern is caused by interference between the epidemics of the two diseases. With this result in mind, we suggest that disease-disease interactions mediated by non-specific immunity (innate immunity) may have a larger effect epidemic patterns than previously suspected.

## DATA, CODE AND MATERIALS

The model is completely reproducible from the description in the main text, and all parameters are made available in the figure captions.

- [1] RE Hope-Simpson and DB Golubev. A new concept of the epidemic process of influenza a virus. Epidemiology and infection , 99(01):5-54, 1987.
- [2] Randy J Nelson and Deborah L Drazen. Melatonin mediates seasonal adjustments in immune function. Reproduction Nutrition Development , 39(3):383-398, 1999.
- [3] S Smith, L Blanton, K Kniss, D Mustaquim, C Steffens, C Reed, A Bramley, B Flannery, AM Fry, LA Grohskopf, et al. Update: Influenza activity-united states. MMWR. Morbidity and mortality weekly report , 64(48):1342, 2015.
- [4] Alicia M Fry, Aaron T Curns, Kathryn Harbour, Lori Hutwagner, Robert C Holman, and Larry J Anderson. Seasonal trends of human parainfluenza viral infections: United states, 1990-2004. Clinical Infectious Diseases , 43(8):1016-1022, 2006.
- [5] Henri Laurichesse, Daniel Dedman, John M Watson, and Maria C Zambon. Epidemiological features of parainfluenza virus infections: laboratory surveillance in england and wales, 1975-1997. European journal of epidemiology , 15(5):475-484, 1999.
- [6] Molly E Counihan, David K Shay, Robert C Holman, Sara A Lowther, and Larry J Anderson. Human parainfluenza virus-associated hospitalizations among children less than five years of age in the united states. The Pediatric infectious disease journal , 20(7):646-653, 2001.
- [7] MMaes, W Stevens, S Scharpe, E Bosmans, F De Meyer, P D'Hondt, D Peeters, P Thompson, P Cosyns, L De Clerck, et al. Seasonal variation in peripheral blood leukocyte subsets and in serum interleukin-6, and soluble interleukin-2 and-6 receptor concentrations in normal volunteers. Experientia , 50(9):821-829, 1994.
- [8] FN Boctor, RA Charmy, and EL Cooper. Seasonal differences in the rhythmicity of human male and female

## COMPETING INTERESTS

We have no competing interests

## AUTHORS' CONTRIBUTIONS

Gorm Gruner Jensen wrote most of the simulations, participated in developing the model; Florian Uekermann helped producing code, and extracted the time-series data shown in figure 6, participated in developing the model; Lone Simonsen conceived the idea of studying the seasonal ordering of infectious diseases; Kim Sneppen coordinated the collaboration, participated in developing the model; All authors participated in designing the study, draft the manuscript, and gave final approval for publication.

## FUNDING

This research has received funding from the European Research Council under the European Union's Seventh Framework Programme (FP/2007 2013)/ERC Grant Agreement n. 740704.

- lymphocyte blastogenic responses. Immunological investigations , 18(6):775-784, 1989.
- [9] Scott F Dowell. Seasonal variation in host susceptibility and cycles of certain infectious diseases. Emerging infectious diseases , 7(3):369, 2001.
- [10] JJ Cannell, R Vieth, JC Umhau, MF Holick, WB Grant, S Madronich, CF Garland, and E Giovannucci. Epidemic influenza and vitamin d. Epidemiology and infection , 134 (06):1129-1140, 2006.
- [11] Pejman Rohani, David J Earn, B Finkenst¨ adt, and Bryan T Grenfell. Population dynamic interference among childhood diseases. Proceedings of the Royal Society of London B: Biological Sciences , 265(1410):20332041, 1998.
- [12] Samit Bhattacharyya, Per H Gesteland, Kent Korgenski, Ottar N Bjørnstad, and Frederick R Adler. Crossimmunity between strains explains the dynamical pattern of paramyxoviruses. Proceedings of the National Academy of Sciences , 112(43):13396-13400, 2015.
- [13] Masashi Kamo and Akira Sasaki. The effect of crossimmunity and seasonal forcing in a multi-strain epidemic model. Physica D: Nonlinear Phenomena , 165(3):228241, 2002.
- [14] Yunxin Huang and Pejman Rohani. The dynamical implications of disease interference: correlations and coexistence. Theoretical population biology , 68(3):205-215, 2005.
- [15] Daniel A Vasco, Helen J Wearing, and Pejman Rohani. Tracking the dynamics of pathogen interactions: modeling ecological and immune-mediated processes in a twopathogen single-host system. Journal of theoretical biology , 245(1):9-25, 2007.

- [16] P Rohani, CJ Green, NB Mantilla-Beniers, and BT Grenfell. Ecological interference between fatal diseases. Nature , 422(6934):885-888, 2003.
- [17] Joaqu´ ın Sanz, Cheng-Yi Xia, Sandro Meloni, and Yamir Moreno. Dynamics of interacting diseases. Physical Review X , 4(4):041005, 2014.
- [18] Anton Camacho, S´ ebastien Ballesteros, Andrea L. Graham, Fabrice Carrat, Oliver Ratmann, and Bernard Cazelles. Explaining rapid reinfections in multiple-wave influenza outbreaks: Tristan da cunha 1971 epidemic as a case study. Proceedings of the Royal Society of London B: Biological Sciences , 278(1725):3635-3643, 2011. ISSN 0962-8452. doi:10.1098/rspb.2011.0300. URL http://rspb.royalsocietypublishing.org/content/ 278/1725/3635 .
- [19] Roberto A Saenz, Michelle Quinlivan, Debra Elton, Shona MacRae, Anthony S Blunden, Jennifer A Mumford, Janet M Daly, Paul Digard, Ann Cullinane, Bryan T Grenfell, et al. Dynamics of influenza virus infection and pathology. Journal of virology , 84(8):39743983, 2010.
- [20] Kasia A Pawelek, Giao T Huynh, Michelle Quinlivan, Ann Cullinane, Libin Rong, and Alan S Perelson. Modeling within-host dynamics of influenza virus infection including immune responses. PLoS Comput Biol , 8(6): e1002588, 2012.
- [21] Rustom Antia and Jacob C Koella. A model of nonspecific immunity. Journal of theoretical biology , 168(2): 141-150, 1994.
- [22] Pengxing Cao, Ada WC Yan, Jane M Heffernan, Stephen Petrie, Robert G Moss, Louise A Carolan, Teagan A Guarnaccia, Anne Kelso, Ian G Barr, Jodie McVernon, et al. Innate immunity and the inter-exposure interval determine the dynamics of secondary influenza virus infection and explain observed viral hierarchies. PLOS Comput Biol , 11(8):e1004334, 2015.
- [23] Neil M Ferguson, Alison P Galvani, and Robin M Bush. Ecological and immunological determinants of influenza evolution. Nature , 422(6930):428-433, 2003.
- [24] Eric Muraille. The unspecific side of acquired immunity against infectious disease: Causes and consequences. Frontiers in Microbiology , 6, 2015.
- [25] Mihai G Netea, Jessica Quintin, and Jos WM van der Meer. Trained immunity: a memory for innate host defense. Cell host &amp; microbe , 9(5):355-361, 2011.
- [26] Mihai G Netea. Training innate immunity: the changing concept of immunological memory in innate host defence. European journal of clinical investigation , 43(8):881-884, 2013.
- [27] Jennifer R Hamilton, David Sachs, Jean K Lim, Ryan A Langlois, Peter Palese, and Nicholas S Heaton. Club cells surviving influenza a virus infection induce temporary nonspecific antiviral immunity. Proceedings of the National Academy of Sciences , 113(14):3861-3866, 2016.
- [28] Vibhuti P Dave, Jane E Allan, Karen S Slobod, F Suzette Smith, Kevin W Ryan, Toru Takimoto, Ultan F Power, Allen Portner, and Julia L Hurwitz. Viral cross-reactivity and antigenic determinants recognized by human parainfluenza virus type 1-specific cytotoxic t-cells. Virology , 199(2):376-383, 1994.
- [29] Joshua B Plotkin, Jonathan Dushoff, and Simon A Levin. Hemagglutinin sequence clusters and the antigenic evolution of influenza a virus. Proceedings of the National Academy of Sciences , 99(9):6263-6268, 2002.
- [30] RM Chanock, RH Parrott, KM Johnson, AZ Kapikian, and JA Am Bell. Myxoviruses: Parainfluenza 1. American Review of Respiratory Disease , 88(3P2):152-166, 1963.
- [31] Caroline Breese Hall, Edward E Walsh, Christine E Long, and Kenneth C Schnabel. Immunity to and frequency of reinfection with respiratory syncytial virus. Journal of Infectious Diseases , 163(4):693-698, 1991.
- [32] Roy M Anderson, Robert M May, and B Anderson. Infectious diseases of humans: dynamics and control , volume 28, chapter 2. Wiley Online Library, 1992.
- [33] Andrei Korobeinikov and Graeme C Wake. Lyapunov functions and global stability for sir, sirs, and sis epidemiological models. Applied Mathematics Letters , 15 (8):955-960, 2002.
Gorm Gruner Jensen Florian Uekermann Kim Sneppen Lone Simonsen
Could short-lasting non-specific immunity explain seasonal spacing of epidemic diseases?
Common respiratory viruses cause seasonal epidemics in sequential patterns.
The underlying mechanisms for this pattern have been debated for some time. For
influenza, contenders include temperature, humidity, and vitamin D levels.
While such seasonal drivers may be sufficient to explain midwinter peaking it
is unclear if other respiratory viruses peaking in spring, summer, or fall have
individual environmental drivers. Here we present a dynamic model of
interacting diseases, in an effort to explain observed seasonal patterns
without requiring multiple seasonal drivers. Our model extends the classical
SIRS-model to include multiple diseases that interact via a short-lasting,
non-specific immunity component. This temporal protection is triggered by
exposure to any epidemic disease and lasts only a couple of weeks. We show that
the inhibiting disease-interaction allows recurrent epidemic behaviour without
any seasonal driving. In the presence of a single seasonal driver our model
parsimoniously predicts epidemic patterns such as sequentially occurring
epidemic diseases. We also present a two-disease simulation reproducing
multiple features observed in time series of parainfluenza strains PIV-3 and
biennial PIV-1 epidemic patterns. Complex epidemic patterns of seasonal
diseases may be explained by non-specific innate or cloub cell (T-cell)
mediated immune responses that result in interaction between unrelated
respiratory epidemic diseases.