V B Greenman
email: j.v.greenman@stir.ac.uk

Pasour

J V Greenman

V B Pasour

Phase control of resonant systems: Interference, chaos and high periodicity

Keywords: eco-epidemiological models, epidemic suppression, subharmonics, noise excitation, dengue fever Control of resonant systems

Much progress has been made in understanding the effect of periodic forcing on epidemiological and ecological systems when that forcing acts on just one part of the system. Much less is known about situations in which several parts of the system are affected. In this case the interaction between the impacts of the different forcing components can lead to reinforcement of system responses or to their interference. This interference phenomenon is significant if some forcing components are anthropogenic for then management might be able to exercise sufficient control to bring about suppression of undesirable aspects of the forcing, for example resonant amplification and the problems this can cause. We set out the algebraic theory when forcing is weak and illustrate by example what can happen when forcing is strong enough to create subharmonics and chaotic states. Phase is the key control variable that can bring about interference, advantageously shift nonlinear response curves and create periodic states out of chaos. The phenomenon in which high period fluctuations appear to be generated by low period forcing is examined and different mechanisms compared in a two-strain epidemiological model. The effect of noise as a source of high period fluctuations is also considered.

Introduction.

Resonance occurs when a disproportionately large response is produced through the excitation of one or more of the natural oscillatory modes of a system by external periodic forcing. Resonance is familiar in engineering systems 3 where it can lead to catastrophic structural failure. Its effect is also apparent in the life sciences where it is thought to be the primary cause of many periodic epidemics in plant, animal and human diseases. These epidemics can have a serious negative impact on the balance of species in wildlife ecologies and on public health. This connection between epidemics and resonance was confirmed in the study of childhood diseases in the 1980s [START_REF] Dietz | Overall population patterns in the transmission cycle of infectious disease agents[END_REF]. Since then resonance has been observed in a wide range of ecological and epidemiological situations [START_REF] Altizer | Seasonality and the dynamics of infectious diseases[END_REF], differing in the period, strength and configuration of the external forcing and in the dynamics of the unforced system.

Seasonality is the most familiar driver of external environmental forcing, affecting both birth and death rates as well as infection transmission and predation strength. Other environmental drivers such as ENSO (i.e. El Nino Southern Oscillation) and rain patterns in Africa and Asia [START_REF] Wichmann | Extinction risk in periodically fluctuating environments[END_REF] have multiannual periods. Environmental forcing is not restricted to such global or regional climate variations. For example, a species, attempting to invade a resident community or web, will see that community or web as the "environment" and will be subject to its often long period fluctuations endogenously or exogenously generated [START_REF] Berryman | The theory and classification of outbreaks[END_REF], [START_REF] Dwyer | The combined effects of pathogens and predators on insect outbreaks[END_REF].

The strength of forcing is also highly relevant in many forcing situations. If forcing is strong enough, subharmonics can be generated whereby the system populations oscillate with a period that is an integer multiple of the forcing period. For example, before vaccination programmes were introduced, measles epidemics occurred typically every two years [START_REF] Dietz | Overall population patterns in the transmission cycle of infectious disease agents[END_REF] even though the forcing period is thought to be seasonal, following the rhythm of the school year [START_REF] Keeling | Seasonally forced disease dynamics explored as switching between attractors[END_REF]. For even stronger forcing the population fluctuations can become chaotic [START_REF] Ireland | The effect of seasonal host birth rates on population dynamics: The importance of resonance[END_REF]. Also of importance is the configuration of the forcing, i.e. which parts of the system are directly affected by the external forcing. For childhood diseases seasonal forcing works predominantly through infection transmission but in other cases other processes can be targeted by the forcing as well. For example for house finches infected with bacterial conjunctivitis (Mycoplasma gallisepticum) breeding occurs in the summer while outbreaks of the infection usually occur in the fall and winter when there is social aggregation [START_REF] Altizer | Seasonality and the dynamics of infectious diseases[END_REF], [START_REF] Hosseini | Seasonality and wildlife disease: how seasonal birth, aggregation and variation in immunity affect dynamics of Mycoplasma gallisepticum in house finches[END_REF]. In contrast, for harbour seals (Phoca vitulina) infected with the phocine distemper virus, breeding and social aggregation occur at the same time, when the seals haul out [START_REF] Altizer | Seasonality and the dynamics of infectious diseases[END_REF], [START_REF] Swinton | Persistence thresholds for phocine distemper virus infection in harbour seal Phoca vitulina metapopulations[END_REF]. When there are multiple components to the forcing they can interact to reinforce responses or, alternatively, to bring about interference between these responses. The spectrum of behaviour from reinforcement to interference is the result of high sensitivity of the amplitude of the population oscillations to the lags between components. There is a half cycle lag in eco-epidemiological systems if maximum predation is in the summer (e.g. the winter spent in hibernation) and maximum prey infection occurs in the winter (due to higher prey density). If instead predation is constant, as it might be in the house finch ecosystem with domestic cats as predators, then there would be a half cycle lag between breeding and infection transmission.

An important case of multiple forcing is when some of the components are the result of human activity. One example of this is the periodic harvesting of a plant or animal population seasonally driven and vulnerable to infection. Such an interaction between anthropogenic and environmental forces can have undesirable consequences if reinforcement rather than interference takes place. Harvesting, for example, can increase the incidence of infection [START_REF] Choisy | Harvesting can increase severity of wildlife disease epidemics[END_REF]. On the other hand, the ability of management to manipulate some forcing components suggests the possibility of controlling the system through "countercyclical" methods, using periodic forcing to offset the harmful effects of the environmental forcing, for example species extinction or large scale disease epidemics. Understanding when the different components of the external forcing will lead to reinforcement and when interference will result is therefore important.

We address this problem in this paper, studying the behaviour of continuous-time models subject to multicomponent periodic external forcing. For weak forcing, when the system can be linearised, the analysis can be entirely algebraic. For stronger forcing we explore through simulation whether well chosen forcing lags are still able to switch the system from reinforcement to interference when subharmonics and chaos are present.

Throughout the discussion the response curve for a system population is used to visualise important aspects of the analysis. To construct this curve the amplitude or the maximum (over a cycle) of the chosen population is plotted against the forcing period p. This is the "natural" way of studying resonance since resonance is primarily concerned with the relationships between three different periods, the forcing period itself, p, the period p 1 of the populations (in response to the forcing) and the (leading) natural period of the unforced system p 0 . This natural period is the period of the damped oscillatory path taken by the (assumed) overcompensating system to bring the system back into equilibrium. This is the dynamical state that is excited by the forcing to bring about resonance. In the simplest cases this excitation shows up as a peak on the response curve when the external and natural periods coincide. Period p 0 is related to and can be calculated from the imaginary part of the leading complex eigenvalue of the Jacobian of the unforced system. (In high dimensional systems there can be more than one natural state that can be excited and hence more than one complex Jacobian eigenvalue and more than one natural period. However the impact of such secondary excitation can be quite small.)

For definiteness we focus on resonant amplification in (eco-)epidemiological systems. The problem is how to suppress the often large scale epidemics that can periodically occur in a population host to a pathogen. To study this problem we use the simplest non-trivial model available, the SI model (S = susceptible, I = infectious population). It has a simple structure without explicit built in delays (e.g. periods of latency or immunity) but is still liable to resonant amplification. Management control will be exercised mostly through host culling.

In the final part of the discussion we focus on one particular feature of the SI model behaviour, the occurrence of seemingly high period fluctuations in system populations when the external forcing has low period. The introduction of a model describing the dynamics of a pathogen with two strains shows that the mechanisms that produce high periodicity for the SI model may also work in more complicated systems. As well, new ways of creating high periodicity are explored, including noise excitation.

Multicomponent forcing.

The fundamentals.

The SI epidemiological model is defined by the equations: while β is the infection transmission constant. For the moment we take P = 0 explaining its significance later.

Equations (1b), (1c) with P = 0 and (1a) redundant define for us Model 1.

To model periodic external forcing we will suppose that there are two components, working through birth rate 'a' and infection transmission 'β'. Precisely:

a = a 0 (1 + δ 1 cos(ωt)), β = β 0 (1 + δ 2 cos(ωt+φ)) (2)
where a 0 , β 0 are average values, δ 1 , δ 2 the forcing strengths, p = 2π/ω the common forcing period and φ the phase between components. The phase can be interpreted as a lag of (-φ/2π), so if φ = -π/2 then β is lagging 'a' by a quarter cycle but if φ = + π/2 then 'a' is lagging β by a quarter cycle. Sinusoidal functions will be used throughout to model periodicity because of the analytical advantages of so doing and because the qualitative features of the dynamics are usually not sensitive to the precise functional form chosen [START_REF] Macdonald | A computational investigation of seasonally forced disease dynamics[END_REF].

With weak forcing (i.e. δ 1 , δ 2 << 1) the model equations can be linearised with explicit formulae given for the amplitudes of the population oscillations (see Appendix B, (B9) and (B10)). These amplitudes are written as a ratio with the denominator (complex) zeros generating the resonance peaks and the numerator moderating or distorting these peaks. If the numerator for one of the populations becomes zero (for some value of the external period p) then the oscillations for that population are eliminated entirely and for all time (for that value of p) leaving the population at its (unforced) equilibrium value. This is an important possibility with significant policy implications.

The conditions for the numerator to be zero, when there are two forcing components, are found algebraically in Appendix B (B7). These are the conditions for the individual responses generated by the two forcing components to have the same magnitude and to be exactly out of phase, reflecting the fact that responses add in a linear system.

These two conditions relate the forcing period p, the phase φ and relative forcing strength θ 0 = δ 1 /δ 2 . Given one of these forcing parameters, the zero amplitude conditions determine the values of the other two if there exists a feasible solution with p positive, which is not always the case.

To illustrate the zero numerator conditions (B7) we analyse Model 1 with the parameters listed as set 1 in we plot the amplitude of this population's oscillations against external forcing period, p). It transpires that there is a zero in the numerator of the infectious population and hence no infectious population oscillations when p = 4.2

(point Z). The suppressive effect of this zero is felt throughout the response curve, in stark contrast to the response curve for phase φ = -π/2 where there is reinforcement rather than interference between the forcing component responses, creating the very sharp resonance peak at natural period p 0 = 2.8. The switch from φ = -π/2 to φ = + π/2 reduces the height of the resonance peak by a factor of 5 even though the numerator zero at Z is well away from the peak position. The effect is even more striking when θ 0 = 2.2 (Fig 1B). Then the positions of the resonance peak and numerator zero, Z, coincide. As a result, the peak maximum (φ = -π/2) becomes a zero minimum (φ = + π/2).

One would not expect that all of the subpopulations of a system have zero amplitude for the same forcing period since this would mean that they all remain at their (unforced) equilibrium values, i.e. there would be no direct evidence of the presence of the external forcing. Algebraically, the zero conditions would be overconstrained and have no solution. For Fig It is not unusual for reliable time series data to exist for only one of the populations from an ecological system of interest. The danger is that from a visual analysis of this one series it might be incorrectly concluded, from the absence of any obvious signs of resonance behaviour, that resonance is not present in the system. As we have noted, it can be suppressed in one population but be present in others.

The time unit used in the analysis need not be "one year". We can rescale to any time unit (with corresponding change in parameter values) in order to show that amplitude zeros can happen for any choice of external period. Our model with the time unit rescaled to "a quarter year" illustrates this fact for seasonal forcing (Fig 1A).

Subharmonics.

With weak forcing we have seen that two (or more) components of the external forcing can destructively interfere and so lead to the suppression of population fluctuations. Does this remain true when the forcing strength is increased sufficiently to trigger subharmonics? Before addressing this question we first review the properties of subharmonics in the simpler case that there is single component forcing on the system. (For a more extended discussion see Greenman et al. [START_REF] Greenman | External forcing of ecological and epidemiological systems: a resonance approach[END_REF].) period. As an example, measles epidemics, before vaccination programmes were started, typically occurred every other year even though the forcing was annual. The interesting question is how subharmonics emerge from a model's response curve as the forcing strength δ is increased. We answer this question in Fig 2 for the SIR measles model, increasing δ from 0.0035 to 0.053 (for details of the model see Appendix C(i)). The qualitative features of this process are in fact common to many epidemiological models, in particular the SI model previously discussed.

We choose the SIR model to construct If the strength of forcing is increased then the response curve becomes distorted (Fig 2B ), the peak P 1 leaning to the right, reflecting "dissipation" in the system [START_REF] Jordan | Nonlinear Ordinary Differential Equations[END_REF]. As a result the simple link between natural and forcing period is broken, the peak no longer occurring when the forcing and natural periods are close in value. The distortion can be so great that the peak overhangs, creating a "breaking wave" profile in analogy with waves breaking on a seashore. This overhang generates two attractors corresponding to large (a) and small (b) fluctuations.

(Points on the dashed curve in Fig 2B indicate the presence of a repeller.) The existence of multiple attractors is of significance, for example, in epidemiology, suggesting that it might be possible to avoid major epidemics by switching to a low amplitude state through a carefully designed control strategy.

Increased forcing not only distorts the original bell shaped peak structure but also introduces new peaks into that curve. The first to appear is, normally, a peak P 2 that emerges from the response curve AP 1 B at the point p = p 0 /2 (Fig 2B ) about which the population oscillates with period, p 1 = 2p, twice that of the forcing period, i.e. this is the first subharmonic. With p 1 = 2p and p ~ p 0 /2 we deduce that p 1 ~ p 0 about the new peak, i.e. the populations are oscillating at or close to the natural period. We will refer to the point, p = p 0 /2, where the peak P 2 first emerges, as the "root" of the peak P 2 and its associated subharmonic. (The root for the basic peak P 1 , reflecting the initial impact of resonance on the dynamics, lies at p = p 0 (Fig 2A).)

As the forcing strength increases, the peak structure about this second peak P 2 , increases in height, in its distortion to the right and in the expansion of the base of the structure (BC in Fig 2D ) about the root (which stays fixed). As a consequence the simple link between forcing and natural periods and between root and peak is lost, as happened with P 1 . Also, the peak structure about P 2 can lean so far to the right that an overhang is formed, thereby creating for itself two attractors (α, β in Fig 2C ).

With yet further increase in the forcing strength more and more subharmonic sections of the response curve appear, each having a "spray" structure (see P 3 , P 4 , P 5 ,…in Fig 2D). Unlike the first subharmonic they are not connected to the baseline response curve ABP 2 CP 1 H 2 D but suspended above it. (The term "spray" extends our analogy of waves breaking on the seashore.) These new subharmonics have period p 1 = np (n integer and greater than two) and "appear" to be rooted on the p axis of the response diagram at points p = p 0 /n.

An additional qualitative change that can occur with strong forcing is the "growth" of period doubling cascades on subharmonic structures. In Fig 2D this phenomenon has already started on the first subharmonic peak structure (about P 2 ). Over its flank BP 2 the population fluctuates with period 2p except on the section (ef) where the fluctuations have period 4p. More period doubling occurs with further increases in forcing, eventually leading to chaos. Period doubling has also started on the 3p subharmonic (P 3 ).

One implication of period doubling is that subharmonics with the same period can be generated in different ways. For example, period 4p oscillations can be obtained after period doubling on the first 2p subharmonic (P 2 )

and also on the 4p subharmonic (P 4 ) before period doubling. One can discriminate between these subharmonics by

Strong forcing with multiple components.

To understand what can happen when there are two or more strong forcing components consider strength is so great that the distortion creates multiple attractor overhangs. In our example not only is the basic peak suppressed but also the emerging 2p subharmonic.

Nonlinear control.

If the two components of forcing are both environmentally driven (e.g. birth rate and infection transmission as in Fig 1) then the forcing parameters δ i , φ, p are all already determined and so there is no guarantee that suppression will take place, i.e. that equations (B7) will be satisfied.

However if one of the components is anthropogenic in origin (e.g. under management control) then there is the possibility of using the freedom to choose forcing parameter values (e.g. To show what can be achieved by such a policy in such a situation consider Model 2 with parameter set 3 of Appendix A with anthropogenic forcing on P and environmental forcing on β: to that achieved at peak points P 2 and P 1 . This variability of structure in a highly distorted response curve suggests an alternative approach to control not based on numerator zeros but rather on the fact that we can change the size and shape of the response curve by altering the phase. With such manipulation we can perhaps arrange that the system is fluctuating not at a peak of the response curve but at one of its low points, between peaks. and the 4p subharmonic 8 times the minimum (the epidemics occurring every 4 years)(point C in Fig 5). Random search of initial conditions suggest that these higher subharmonics can be as easily accessed as the low amplitude seasonal oscillations (point A in Fig 5). Noise can switch the system between these 3 possible modes as has been found in childhood diseases [START_REF] Greenman | External forcing of ecological and epidemiological systems: a resonance approach[END_REF]. So in this situation the change in phase (2.4 ->(-0.8)) would not be recommended. However a change to a different less optimal but still advantageous phase can avoid this difficulty.

P = P 0 (1 + δ 1 cos((2πt/p)+φ)), β = β 0 (1 + δ 2 cos(2πt/p)) (3

A realisation of this idea is given in

The natural period of the model of Fig 5 (p 0 = 2.5) warns us of the danger from these subharmonics. Its value indicates that the first (2p) subharmonic is not likely to be accessible with seasonal forcing (since root p 0 /2 > 1) but higher subharmonics (3p, 4p) could be if the forcing is strong enough (since their roots p 0 /3, p 0 /4 are both less than one).

3. Multiple forcing with chaotic solutions. interval in which the behaviour appears to be periodic with period 15. During these cycles the susceptible population rises from a low level after a major epidemic up to its carrying capacity, K = 2.5, while the number of those infected is so small that they have no effect on this recovery process. The pathogen needs this amount of time to re-establish itself. Then its rapid exponential growth overwhelms the sizeable susceptible population which collapses in the face of overwhelming odds. The consequent scarcity of susceptibles available to the pathogen means that the infectious population collapses as well. The period of 15 is roughly the time the susceptible population needs to recover from the previous epidemic.

Over a much longer time interval the system is chaotic as evidenced by its trajectory plotted in state (S, I) space in Fig 7B . So it is only over short intervals that the time series appears periodic. To study further the properties of this chaotic solution we transform the problem from a continuous into a discrete time problem by taking (stroboscopic) Poincare sections of the state space trajectory after every p units of time [START_REF] Strogatz | Nonlinear Dynamics and Chaos[END_REF]. These points form a discrete time trajectory in state space from which we can extract the Fourier power spectrum to confirm, in the first instance, whether or not the dynamics is chaotic [START_REF] Seydel | Practical Bifurcation and Stability Analysis[END_REF]. If the dynamics is periodic or quasiperiodic there will be a single or multiple spikes defining the left hand part of the spectrum but if chaotic then there will typically be a "sea of noise" from which emerges power spikes, suggesting chaotic behaviour about an underlying (quasi)periodic structure. The position of the dominant spike will provide us with a measure of the periodicity of the chaotic timeseries. Consider, for example, the Fourier power spectrum for p = 1 (Fig 7C) which shows chaos. Its left hand dominant peak lies at position 340 indicating a discrete time "chaotic period" of 5000/340 ~ 15 and hence a "chaotic period" of 15p = 15 (with p = 1) in continuous time. This is also the period of the short term "periodic" time-series shown in Fig 7A . (Note that the figure of 5000 quoted in the chaotic period calculation is the number of points on the system trajectory generated by Poincare section.)

Of importance are two indices to characterise the chaotic behaviour, the maximum infectious population reached over the chaotic trajectory (to be referred to as the "chaotic maximum" (e.g. sufficiently long period of time to enable the chaotic maximum or the window maximum to be identified. The periodicity (chaotic or otherwise) can be determined from the Fourier power spectrum. The resolution in p is adjusted as necessary to capture the detailed structure.

Chaotic period and chaotic maximum.

There are two important properties of the chaotic part AB of the response curve for φ = π in Fig 6 that should be noted. First, the "chaotic period" in continuous time is roughly constant over AB as can be checked by simulation. Not only is the chaotic period near constant but it is near constant at a high value. For example, for p = 1 the chaotic period is 15 times as large as the forcing period. As we previously observed this high value reflects the recovery time from a major epidemic. It is not related to the natural period of the system, i.e. the imaginary part of the leading complex eigenvalue which, in this example, has value 1.86, but to its real part which measures the decay (recovery) time back to equilibrium.

The second important feature of this chaotic behaviour is that the chaotic maximum remains high throughout the chaotic regime, always above 1.25 for 1 < p < 9 (Fig 6 ) whereas the unforced point equilibrium has value 0.008.

As well, there is a slow upward drift in the chaotic maximum as p is increased, rising from 1.25 to 2.2 at p = 8. This is due to the progressive overshooting of the carrying capacity, evident in the comparison between the susceptible time series in Figs 7A,D. Whether this chaotic maximum can be reduced by phase change is a question we answer in the next section.

As a point of interest, the behaviour we see in Fig 7A , of a rise to equilibrium of the resident population that enables "invasion" by the pathogen which then leads to collapse of both resident and invader populations to restart the cycle, has some similarities to the dynamics of the (unforced) prey-parasitoid-hyperparasitoid system studied by Ives and Jansen [START_REF] Ives | Complex Dynamics in Stochastic Tritrophic Models[END_REF]. Approaching equilibrium for the prey-parasitoid subsystem enables the hyperparasitoid to invade, causing the resident subsystem to respond and exclude the invader, leading to a repeat of the cycle. For this tritrophic system the period of the cycle has the exceptionally high value of 250 years.

The effect of phase change.

To determine by how much the chaotic fluctuations can be reduced by phase change we also plot in Fig 6 the response curves for φ = -π/2 (top (dashed) curve) and φ = + π/2 (bottom (dashed) curve). As with the middle curve (for φ = π), the top curve is a continuum of chaotic states punctuated by periodic windows of varying width and depth as the system finds periodic solutions with often much lower amplitude fluctuations. Beyond p = 9 the system is periodic and the maximum is reached around p = 10. The bottom curve, for φ = + π/2, also has a sequence of periodic windows within an otherwise chaotic continuum but the windows are more modest in width and depth. The states are periodic beyond p = 11 with the amplitude still rising.

These two curves show how much can be gained if the phase is switched from φ = -π/2 to + π/2, e.g. when birth rate or infection transmission is managed. For example for p = 6 the chaotic maximum can be reduced by a factor of 4 by the phase shift while for p = 1 the factor is 1.8 with a bonus that the switched state (φ = + π/2, p = 1) is periodic (with period 13), i.e the behaviour becomes predictable. In general, if circumstances are such that achieving a periodic (rather than chaotic) state is advantageous then it might be possible to find a phase φ where the reduction in fluctuation is less but the switched state is periodic.

High period fluctuations.

As one outcome of our analysis we have identified various ways in which the populations of a forced system can have a much higher period than that of the external forcing. This is the topic that we will now focus on for the rest of the discussion.

Subharmonics combined with periodic doubling was the first mechanism for generating high periodicity that we studied. However to bring this about the forcing, in many cases, needs to be so strong that it triggers chaos first but, as the SI model shows, chaos can itself be the source of high period output. This can happen within the windows of periodicity that punctuate an otherwise chaotic continuum as is the case, for example, when input period An interesting question is whether in more complex models these mechanisms can still be present. One might expect a positive answer since the subharmonic phenomenon and the occurrence of periodicity within chaos are not rare occurrences in nonlinear ecological and epidemiological systems. To check out this supposition we now study the dynamics of a model much larger than the SI model where a single population is host to two strains of a pathogen with only partial protection against secondary infection. A study of this model will also reveal new high period mechanisms and will connect our discussion of the high period phenomenon with a specific problem in epidemiology.

4.

High periodicity in a two strain epidemiological model.

Model solutions.

dominant strains of dengue fever. This second application is of particular interest since field data suggests that epidemics in each strain occur roughly every 10 years even though the external driver is thought to be seasonal.

The model has nine states defined by whether an individual is in an S; I; R (susceptible; infectious; recovered) state with respect to each strain [START_REF] Kamo | The effect of cross-immunity and seasonal forcing in a multi-strain epidemic model[END_REF]. Infection transmission is density dependent and seasonal environmental forcing acts through the model's transmission constants.

The key parameter of the model, σ, measures the degree of protection from secondary infection.

To simplify the model, symmetry between strains is assumed, reducing the model to the five equation set listed in Appendix C(iii) where full details of the model can be found.

Simulation of the model shows that chaos sets in for surprisingly low values of the forcing strength, suggesting that high periodicity within chaotic dynamics is likely. This is confirmed in higher with epidemics in each strain every 35 years. However these solutions are unstable and so are not sustained over the long term. Even so these solutions can still be of interest since it may be that the field data also reflects transient behaviour. It is in fact possible to generate a stable high period subharmonic solution without noise but this requires higher infection transmission (see parameter set 8b). In Fig 8D the forcing is seasonal and it is the 4p subharmonic that is involved, as can be established by detailed numerical examination. For the forcing strength chosen and for p =1 this subharmonic has period doubled, generating period 8 oscillations for the susceptible population. However the two strains being out of phase means that the period is further doubled for the infectious populations, yielding the period 16 oscillations of Fig 8D, with several subepidemics present.

Noise excitation.

The dengue field data suggests epidemics of period 10. The dengue model parameters [START_REF] Adams | Cross-protective immunity can account for the alternating epidemic pattern of dengue virus serotypes circulating in Bangkok[END_REF] indicate a leading natural period of 10. There would be a simple explanation for this "coincidence" if the forcing also had periodicity of around 10. We would be back to straightforward resonance excitation. But the forcing is thought to be seasonal. However there is a possible way around this by thinking further about the impact of noise on a dynamical system. Since (white) noise is a homogeneous combination of all forcing frequencies its presence should, in principle, excite all the natural modes of a system, in particular the leading mode. This phenomenon is evident in Fig 9 where we plot a typical time-series for the dengue model with noise (rather than seasonality) being the dominant external force. In this figure we observe noise induced bursts of fluctuations with the period between epidemics "locked" into the leading natural period of 10. It is the amplitude rather than the period that mostly captures the randomness of the forcing. In other time-series realisations there is frequent switching between the two natural modes with periods 6 and 10.

Noise can also be applied to simpler single strain epidemiological models subject to seasonality but, certainly for childhood diseases (such as measles and rubella) and the examples of the SI model we have considered, the natural period is not too distant from p = 1, the seasonal forcing period. So the noise induced generation of high periods (compared to p = 1) that we have been discussing for dengue will not work for these single strain models even when noise is the dominant external force because there are no high value natural periods to excite. However in trophic (predator-prey) models the leading natural period typically takes a high value and noise excitation can be successful.

Analysing noise.

The ability of noise to excite a system's natural modes can be demonstrated theoretically if the model is in discrete-time and if the forcing is weak enough for the model to be linearised. Then there exists an explicit formula (B9) (Appendix B), derived from standard statistical theory [START_REF] Jenkins | Spectral Analysis and its Applications[END_REF], [START_REF] Reinsel | Elements of multivariate time series analysis[END_REF], giving the frequency (power) spectrum of the population time-series. Typically the largest contributions come from the frequencies (periods) corresponding to the natural modes and these modes define the basic structure of a randomly generated time-series. A good example of this is the dynamics of North Atlantic cod modelled by Bjornstad et al. [START_REF] Bjornstad | Cycles and trends in cod populations[END_REF]. The data exhibits both short and long term fluctuations, a feature that has been explained as the simultaneous excitation by noise of the two natural modes of the model [START_REF] Greenman | The frequency spectrum of structured discrete time population models: its properties and their ecological implications[END_REF].

The theory of noise excitation of biological systems in continuous time is far more complex [START_REF] Gard | Introduction to stochastic differential equations[END_REF], [START_REF] Cai | Stochastic analysis of the Lotka-Volterra model for ecosystems[END_REF] if one assumes uncorrelated random variation. A perhaps biologically more realistic assumption is to follow Adams et al. [START_REF] Adams | Cross-protective immunity can account for the alternating epidemic pattern of dengue virus serotypes circulating in Bangkok[END_REF] and suppose strong short term correlation which can be approximately modelled by assuming forcing that changes (randomly) only after a given time interval. In this way one can model seasonal or annual variation rather than inconsequential second by second random changes (Appendix C(iii)). This approach allows standard differential equation simulation to be used (Fig 9).

Verification.

Noise excitation of high period epidemics in a model with high natural period adds to the range of possible mechanisms that should be considered in analysing high period phenomena in particular ecological and epidemiological systems. As with the other mechanisms discussed, the focus has been on explaining the periodicity rather than other features of the field data such as the magnitude of the epidemics and the existence and pattern of any subepidemics. These further tests would help in discriminating between the various possible explanations as would, as always, more field data.

Discussion.

Several messages emerge from our analysis of externally forced systems prone to resonant amplification. First is the importance of phase relationships when there are multiple forcing components. We found in our models that varying the phase can switch the response for one of its populations from reinforcement to interference, achievable both when the forcing is weak (Fig 1 ) and when the forcing is strong enough to trigger subharmonics (Fig 3) and chaotic states (Fig 6). The ability to switch between reinforcement and suppression is of limited interest if the forcing is entirely environmental (since the phases are then pre-determined) but can be of practical interest if one or more of the forcing components is under management control. In the latter case the possibility exists of significantly reducing or even eliminating the fluctuations in at least one of the system's populations, of advantage if these fluctuations are undesirable or costly. The application we have studied in detail is the suppression of epidemics (i.e. fluctuations in the infectious populations) in order to reduce the infection to much lower levels. Also of relevance is the problem of large fluctuations posing a threat to a population's survival by driving the population to dangerously low levels. The risk for management is that if the system dynamics is not sufficiently understood the phase controls might be set at values that would worsen rather than improve the situation.

The second observation is that extra caution has to be shown when looking for resonance in a limited amount of data. If, for example, the only reliable data at hand is for just one of the populations then the lack of evidence of resonance in that data does not necessarily mean that resonance is absent from the system. It might be suppressed in the population for which there is data but not in other populations. In principle one should be able to reconstruct the hidden structure if sufficiently accurate and sufficiently long time series data is available for just one of the populations [START_REF] Takens | Detecting strange attractors in turbulence[END_REF] but this is seldom the case.

An important issue emerged from the discussion on reinforcement and suppression of chaotic fluctuations, namely how it is possible that forcing with a low period (e.g. seasonality) can generate a seemingly high period output in those situations where the forcing is unable to trigger the required subharmonic. This finding brings together a range of issues. The dynamics might appear to be periodic in the short term but this behaviour could be transient or short lived. On the other hand long term periodicity could be transient if it were not sustained by the presence of noise. We also noted that noise alone can excite the natural oscillations of a system, however high the natural period might be while periodic forcing might not be able to do so. These observations indicate the importance of analysing, as well, the unstable and short-lived behaviour of models in explaining the field data when it is not clear whether that data reflects transient or long term behaviour [START_REF] Hastings | Transient dynamics and persistence of ecological systems[END_REF] There is also the possibility that even if there is some periodic or stochastic forcing present the primary driver of population fluctuations is endogenous, i.e. self-excitation takes place [START_REF] Koelle | Disentangling Extrinsic from Intrinsic Factors in Disease Dynamics: A Nonlinear Time Series Approach with an Application to Cholera[END_REF]. An example of this occurs in the modelling of dengue. Recker et al. [START_REF] Recker | Immunological serotype interactions and their effect on the epidemiological pattern of dengue[END_REF] constructed a model that takes into account antibody-dependent enhancement, increasing both susceptibility and infection transmission during secondary infection. In this model high period fluctuations can happen without external forcing. Increasing the enhancement parameters leads to system instability. Whether dengue epidemics are driven exogenously [START_REF] Adams | Cross-protective immunity can account for the alternating epidemic pattern of dengue virus serotypes circulating in Bangkok[END_REF] or endogenously [START_REF] Recker | Immunological serotype interactions and their effect on the epidemiological pattern of dengue[END_REF] remains a matter for further discussion.

We have based much of the discussion on the simple SI model but investigation of models with more structure (e.g. the SEIR epidemiological model with latency (E) and recovery (R)) shows similar behaviour. The importance of additional structure is that new opportunities arise for dynamic system control. This is illustrated by previous work on structured discrete time ecological systems [START_REF] Greenman | The impact of environmental fluctuations on structured discrete time population models: Resonance, synchrony and threshold behaviour[END_REF] which shows that suppression can in fact occur when there is only one forcing component, the interference happening internally between the "energy" flows created by the forcing on different pathways through the system. The best way to study this mechanism is to work with the network graph, representing the structure of the system and its external connections under weak forcing. This technique was used to illustrate a second way in which fluctuations can be suppressed when there is single component forcing and unique pathways. Suppression can occur through "bottlenecks" on pathways and through "absorption" into subnetworks of the graph. These mechanisms suggest further ways in which management can intervene to achieve objectives. Combining control through external forcing and through intervention in the internal structure of the unforced system can open up a wider range of control options.

We have discussed countercyclical control methods that can be employed to suppress disease epidemics in an externally forced eco-epidemiological system. In some situations, for instance in agriculture, this may be sufficient since removing the disease completely may not be cost effective. The major economic damage is usually caused by sudden epidemics of a plant disease devastating a seasonal crop. When it is deemed essential to exclude a pathogen completely, then using countercyclical control measures might be one way to help achieve this.

That this is a feasible approach has been shown in a range of realistic examples where substantial reduction in the effort involved in excluding a pathogen can be made [START_REF] Greenman | The impact of environmental fluctuations on structured discrete time population models: Resonance, synchrony and threshold behaviour[END_REF], [START_REF] Greenman | Environmental forcing, invasion and control of ecological and epidemiological systems[END_REF], [START_REF] Greenman | Pathogen exclusion from eco-epidemiological systems[END_REF]. These examples show that the concepts of reinforcement and interference and the use of phase as a control variable are just as important in solving the exclusion problem as they are the epidemic suppression problem. We hope to add to these numerical studies by using an analytic approach based on the rare invader approximation [START_REF] Greenman | Pathogen exclusion from eco-epidemiological systems[END_REF] to obtain general results for the exclusion of a species. In other recent work in this area [START_REF] Bacaer | Resonance of the epidemic threshold in a periodic environment[END_REF], [START_REF] Bacaer | Periodic matrix population models: growth rate, basic reproductive number and entropy[END_REF] much of the emphasis has been on mathematical rigour, constructing an R 0 threshold index under forcing for example, rather on the exclusion control problem in ecological models experiencing multi-component forcing. Set 3: (a 0 , b, α, γ; K, β 0 ; P 0 , c 1 , c 2 ; δ 1 , δ 2 ) =

(1.5, 0.5, 73.0, 0.5; 1.5, 79.7; 0.2, 1.0, 1.0; 0.03, 0.51).

Set 4: (a 0 , b, α, γ; K, β 0 ; P 0 , c 1 , c 2 ; δ 1 , δ 2 ) = (1.5, 0.5, 73.0, 0.5; 1.375, 79.7; 0.2, 1.0, 1.0; 0.05, 0.50).

Set 5: (a 0 , b, α, γ; K, β 0 ; P 0 , c 1 , c 2 ; δ 1 , δ 2 ) = (1.5, 0.5, 73.0, 0.5; 1.5, 79.7; 0.2, 1.0, 1.0; 0.05, 0.50).

Set 6: (a 0 , b, α, γ; K, β 0 ; P 0 , c 1 , c 2 ; δ 1 , δ 2 , φ) = (1.5, 0.5, 10.9, 0.0; 4.5, 12.0; 0.5, 0.5, 2.5; 1.0, 1.0, 7π/4).

Set 7: (a 0 , b, α, γ; K, β 0 ; P, δ 1 , δ 2 ) = (1.0, 0.5, 35.8, 0.0; 2.5, 39.1; 0.0, 1.0, 1.0).

Set 8a: (β 0 , γ, μ, σ; p, δ) = (150, 50, 0.01, 0.5; 1, 0.1).

Set 8b: (β 0 , γ, μ, σ; p, δ) = (350, 50, 0.01, 0.5; 1, 0.1).

For example for Model 1, described by equations (1c), (1b) and ( 2 where m ij is the (i,j) element of the Jacobian matrix: Note that in forcing the birth rate 'a' in Model 1 we assume that s = (r/K) = (a-b)/K in (1a), (1c) does not vary. Put another way, both r and K both vary with the seasons and we suppose, for simplicity, that their variations are in phase and cancel out [START_REF] Bowers | Evolutionary branching/speciation: contrasting results with explicit or emergent carrying capacities[END_REF]. Finally note that by making equation (1a) redundant the forcing components can be associated with different equations since 'a' appears only in (1c) and β only in (1b). This has advantages in calculating amplitude zeros.

M = r 0 (1-2H 0 /K) -α β 0 I 0 β 0 S 0 -d -β 0 I 0 ⎡ ⎣ ⎢ ⎢ ⎤ ⎦ ⎥ ⎥ (B2).
(ii) Finding amplitude zeros.

From (B1) with i = √(-1), I the identity matrix, x = Re[z] and z = z 0 e iωt we have:

z 0 = Σ δ l exp(iφ l ) (iωI-M) -1 b l (l = 1, …., m) (B4)
using the (complex) exponential form for trigonometric functions. With B(ω) := (iωI-M) -1 , the jth population complex amplitude is given by: z 0j = ΣΣ A ljk (ω)exp(i(φ l + ψ jk ))

hence:

x j = ΣΣ A ljk (ω)cos(ωt + φ l + ψ jk ) (B5). Chaotic behaviour for Model 1 with forcing on 'a' and 'β' and phases ± π/2, π (parameter set 7).

Subharmonic periodic windows with periods as stated. (q indicates a subharmonic in a higher family of subharmonics (cf Q 1 , Q 2 in Fig 2D ) while ch = chaos.) High period fluctuations for the infectious populations (y i ) with strain i in the two strain model.

  dS/dt = a H -s H 2 -b S -β S I + γ I -c 1 P S (1a) dI/dt = β S I -d I -c 2 P I ( 1 b ) dH/dt = r H -s H 2 -α I -c 1 P S -c 2 P I (1c). S; I; H = S + I denotes the susceptible; infectious; total population while a; b; r = a-b measures per capita birth rate; mortality; net growth rate at low population levels. Self-regulation is introduced through the carrying capacity K, where s = r/K. Parameter α; γ; d = b + α + γ describes virulence; recovery rate; loss rate from the infectious state

  Appendix A. In Fig 1A with θ 0 = 1.5 and φ = + π/2, we plot the response curve for the infectious population, I, (i.e.

  Fig 2 since it is possible to find in this model all the structural features needed in the rest of the discussion in just one response curve (δ = 0.053, Fig 2D). (Note that we plot the maximum reached by the infectious population against the period p. Amplitude is an inappropriate measure when nonlinearities distort the population oscillations.) If the forcing is weak then the response curve has an almost symmetric bell shape centred about the system's natural period p 0 indicating resonant amplification when the natural and external periods match and excitation of the natural mode of the unforced system takes place (AP 1 B in Fig 2A). The populations oscillate in response to the forcing with the same period p, as indicated by the label "1p" on the response curve in Fig 2A.

  Fig 2D also shows the emergence of two other types of dynamic behaviour. First, there is a new family of subharmonics appearing, adding to the ways in which a given period can be generated. The Q 1 "spray" subharmonic

Fig

  Fig 3A where we show the infectious population response curve for the model of Fig 1B with forcing strength at δ 1 = 0.11, δ 2 = 0.05 (θ 0 = 2.2) and phases chosen for comparison to be ± 1.4. (The full parameter set is given as set 2 in Appendix A). When φ = -1.4 the two forcing

  strength and lag) to achieve management objectives (e.g fluctuation suppression). To study this control problem consider Model 1 with P now taken to be non-zero and interpret it as culling of (i.e. predation on) susceptible and infectious prey, thereby creating for us Model 2. One application is "harvesting" of game birds infected with a parasite and environmentally forced through seasonality in infection transmission. Conservation management has the task of deciding when harvesting should be allowed and at what level. If the primary objective is to reduce infection, then management might consider the use of countercyclical control through culling to counter the fluctuations in infection transmission.

  Fig 3B illustrates the fact that the reduction that might be achieved in the fluctuations of

2 . 5

 25 Fig 4A where we plot part of the infectious population response curve for Model 2 (and set 4 parameter values (Appendix A)) when, as before, forcing is only on infection transmission β and predation P (with common period p). This figure shows the first subharmonic peak when the predation lags infection transmission by a quarter cycle (φ = -π/2) and when it leads by a quarter cycle (φ = π/2). Increasing the lag (-π/2 -> + π/2) shrinks the upper response curve A 1 P 2-, in particular narrowing its base which straddles its (phase independent) root at p = 1.05. (By base we mean the interval BC in Fig 2D on which the first subharmonic P 2 "stands".) This phase manipulation results in subharmonic A 2 P 2+ and a reduction in the fluctuations by 67% (from F -> G) when the forcing is seasonal. The reduction in Fig 4B (based on parameter set 5) is even greater. If we again reduce the phase from a lag to a lead of a quarter cycle the fluctuations are reduced by 77% for p = 1.05 (from F -> G). The difference in this second case is that we are at a low point (G) to the right of the subharmonic peak P 2+ in Fig 4B but to the left in Fig 4A. For p = 1.0 in Fig 4B we have multiple attractors so shocking the system (e.g. by appropriately manipulating the initial susceptible and infectious populations) can take us to the lower solution at B. Higher subharmonics. One possibility we did not mention in Fig 4A is that, although the first (2p) subharmonic is not accessible when φ = π/2, higher subharmonics (e.g. 3p, 4p) might be if the forcing is strong enough (Fig 2D). For the model of Fig 4A this turns out not to be the case but it is the case for Model 2, with parameter set 6, with seasonal forcing on birth rate and predation and with forcing strengths at the maximum value of one. For this model, simulation shows that for p = 1 the infectious population fluctuations are minimised when φ = ( -0.8) and maximised when φ = 2.4, the fluctuations being seasonal in line with the forcing in both cases. However at the value of φ for which the fluctuations are minimised there are multiple attractors (Fig 5) since the 3p and 4p subharmonics can also be accessed. The 3p subharmonic has fluctuations 4.7 times those at the minimum (but the epidemics occur every three years rather than annually)(point B in Fig 5)

3. 1

 1 The structure of chaotic solutions.Dramatic changes in the magnitude of population fluctuations can be brought about by change of phase not only when populations oscillate periodically (Figs 1, 3, 4) but also when the dynamics are chaotic. To show this, consider again Model 1 (with no predation) forced through birth rate 'a' and infection transmission β (see[START_REF] Altizer | Seasonality and the dynamics of infectious diseases[END_REF]) with virulence and transmission scaled high enough to trigger chaos (parameter set 7, Appendix A). With φ = π (i.e. a half cycle lag), the response curve for this system is shown as the middle curve in Fig 6.This reveals a near continuum of chaotic states, AB, from p = 1.0 to 9.0, interspersed with windows of periodicity[START_REF] Nagashima | Introduction to Chaos[END_REF] of varying width and depth. For p below 1.0 the chaos becomes unstable, in the long term collapsing to low amplitude oscillations with population and forcing periods equal. At p = 9.0 chaos seamlessly merges into periodic oscillations with period 2p.The windows of periodicity have an interesting structure. There is an (open) interval about the window minimum where the period is as specified. As one "climbs" away from the minimum, period doubling takes place leading to chaos near the "rim" of the window. The windows are therefore an integral part of the chaotic structure and not a separate component of the dynamics. (A search was carried out for additional but separate attractors but none were found for this SI model but are no doubt present for more complex models.)The source of the chaos in Fig6is the period doubling that takes place on each subharmonic peak structure as forcing strength is progressively increased. The start of this process for the first subharmonic (P 2 ) is shown in Fig2Das section (ef) on its flank BP 2 . Along (ef) the population period is 4p rather than 2p. Further increase in forcing expands section (ef) to accommodate a new section within (ef) where the population period is now 8p through more period doubling. This process continues until chaos is generated within a much expanded (ef).The intervals of chaos so created on this and other peak structures eventually merge, with increased forcing, to form the near continuum of chaotic states (Fig6). Being chaotic does not mean that the response curve of Fig 6 (with φ = π) is without order. To understand the dynamics of what is happening we show in Fig 7A the time-series for populations S, I for p = 1.0 over a time

  Fig 7B)) and the "chaotic period" calculated from the Fourier power spectrum (e.g. Fig 7C) and compared with the near periodic time series (e.g. Fig 7A). The basic procedure in constructing the response curve AB in Fig 6 is to generate the (S, I) trajectory over a

  This constancy is illustrated by a second example, for p = 7.4, shown in Figs7D,E,F. This chaotic solution (Fig 7E,F) has similar properties to the first example with p = 1, including intervals of near periodicity (Fig 7D). The chaotic period is 2500/1140 = 2.2 in discrete time (Fig 7F) and 2.2p ~ 16 in continuous time, consistent with the period of the associated time series in Fig 7D. The periodic windows within the chaotic region AB in Fig 6 also contribute to this near constancy. As examples, the 3p window is centred at p = 6 generating period 18 = 6x3 and the 6p window at p = 3 generating period 18 = 3x6, in continuous time.

is 1 and output period 13 with phase π/ 2 (

 2 Fig 6). Also, short intervals of near periodicity (with high period) can occur within an otherwise chaotic trajectory (e.g. Figs 7A,D), so a limited time series might suggest periodicity although the dynamics is fundamentally chaotic. Further, not only can the near periodicity be short lived but also the survival of the chaotic trajectory itself, i.e. the trajectory is in fact created by a weak repeller which acts like a weak attractor in the short term. An example of this doubly transient behaviour was found in Fig 6 when φ = π and p is below 1. In this case large amplitude chaotic fluctuations eventually collapse into low amplitude low period oscillations.

Figs

  Figs 8A, 8B which show solutions of the model that have brief intervals of near periodic

Fig

  Fig 8C shows another solution (for parameter set 8a) with each strain exhibiting 10 year

  more clearly the other dynamical modes. Variations in γ, K and especially forcing strengths were used to explore the SI resonance region and identify the control mechanisms discussed.Further research is required to produce a more coherent picture of the dynamics. Listed below are the parameter sets used. (set 1; (set 2; set 3); (set 4; set 5); set 6; set 7; set 8a, b) relates to (Fig 1; (Fig 3A; Fig 3B); (Fig 4A; Fig 4B); Fig 5; Fig 6; Fig8). Set 1: (a 0 , b, α, γ, K, β 0 , P) = (1.0, 0.5, 18.25, 0.0, 2.0, 19.92, 0.0). Set 2: (a 0 , b, α, γ, K, β 0 ; δ 1 , δ 2 , P) = (1.0, 0.5, 18.25, 0.0, 2.0, 19.92; 0.11, 0.05, 0.0).

  ) with P = 0 and m = n = 2, the linearised equations are: dx 1 /dt = m 11 x 1 + m 12 x 2 + a 0 H 0 δ 1 cos(ωt) d x 2 /dt = m 21 x 1 + m 22 x 2 + S 0 I 0 β 0 δ 2 cos(ωt + φ)

H 0 ,

 0 I 0 , S 0 are the population equilibrium values and x = [x 1 , x 2 ] = [H -H 0 , I -I 0 ] T (T denotes the transpose operation). Also r 0 = a 0 -b. The linearised equations are of the form (B1) with: b 1 = [a 0 H 0 , 0] T , b 2 = [0, β 0 S 0 I 0 ] T ( B 3 ) .

  Fig 9, h = 0.5 and a = 0.25.

Figure 1 :

 1 Figure 1: Infectious population response curves for weak forcing.

Figure 2 :

 2 Figure 2: Response curves with subharmonics. The evolution of the SIR model response curve (with forcing on β) when (A) δ = 0.0035 (B) δ = 0.0095 (C) δ = 0.0125 (D) δ = 0.053. Pairs (a, b), (α, β) indicate multiple attractors. Dashed

Figure 3 :

 3 Figure 3: Interference when system is nonlinear.

Figure 4 :

 4 Figure 4: Shifting subharmonic peaks with phase change.

Figure 5 :

 5 Figure 5: Higher subharmonics.

Figure 6 :

 6 Figure 6: Response curves for chaotic fluctuations.

Figure 7 :

 7 Figure 7: Near periodicity within chaos. Aspects of chaotic behaviour for the model of Fig 6 with φ = π for p = 1 and p = 7.4. (A, D) short

Figure 8 :

 8 Figure 8: High period output generation.

(Figure 9 :

 9 Figure 9: High period output generation

  Figure 1

  Figure 2

  Figure 3

  Figure 4

x 2 = A 121 (ω)cos(ωt + ψ 21 ) + A 222 (ω)cos(ωt + φ 2 + ψ 22 )

Acknowledgement.

We would like to thank the University of Stirling, through the good offices of Leslie Smith, for the funding, direct and indirect, that made this research possible.

Mathematical Appendices:

(i) Models 1 and 2:

To find a point of entry into the resonance region of the SI model we started with a model parametrised for Fox Rabies [START_REF] Anderson | Population dynamics of fox rabies in Europe[END_REF], [START_REF] Ireland | The effect of seasonal host birth rates on population dynamics: The importance of resonance[END_REF]. This model has the disadvantage that virulence α and infection transmission β take "extreme" values and lead to chaos overwhelming the dynamics (unless this tendency is offset by predation). So, where necessary, α, β were scaled down to see Suppose after linearization about a stable point equilibrium the model equations take the form:

where the summation is over index l taking values 1, 2, …, m while M is the nxn system Jacobian matrix evaluated at the equilibrium. Vector x is the vector of population deviations from their equilibrium values. The constant vector b l in the forcing term of (B1) has zero elements for those equations that do not contain the lth forced parameter. Its non-zero elements can be read off from the form of the linearised equations. (B6)

where φ 2 = φ if the A ljk terms are both positive or φ 2 = φ + π otherwise, to absorb the negative sign. From (B6) it is clear that x 2 will have zero amplitude if the magnitudes of the individual terms in (B6) are equal, i.e. A 121 = A 222 or more clearly:

and if they are exactly out of phase:

So to apply these zero amplitude conditions to Model 1 we need to calculate the magnitudes and arguments of the complex numbers B 21 , B 22 . Note that for x 1 and x 2 to have simultaneous zeros we need to satisfy four equations for the three forcing parameters θ 0 , φ, p. In general this over constrained problem will have no solution unless we consider an internal parameter to be a fourth control variable.

(iii) The cofactor amplitude formula in the general case. with summations over k and m. Here:

with summations over l and n and '*' indicating complex conjugation.

Formula (B9) also holds for stochastic forcing if we interpret (B9) as the "power" and ω the frequency in the frequency spectrum for the jth population and S km in (B10) as the covariance matrix for the random variables acting on the different populations. The (complex) zeros of the denominator Δ(ω) in (B9) generate the peaks in the power (frequency) spectrum that denote excitation of the natural modes of the unforced system by stochastic forcing [START_REF] Jenkins | Spectral Analysis and its Applications[END_REF], [START_REF] Greenman | The impact of environmental fluctuations on structured discrete time population models: Resonance, synchrony and threshold behaviour[END_REF].

(ii) Subharmonics in Fig 2D:

A subharmonic can be characterised by two positive integers (m, n) so that the period can be expressed as np and the root lies at (m/n)p 0 . The number of subepidemics over a cycle is (m -1). For the first family of subharmonics (P i ), index m equals 1 and there are no subepidemics. In contrast, Q 2 in Fig 2D has n = 6, m = 3 (i.e 2 subepidemics) while Q 1 has n = 3, m = 2 with 1 subepidemic. (See [START_REF] Greenman | External forcing of ecological and epidemiological systems: a resonance approach[END_REF] for further details.)

(iii) The dengue model [19]:

The nine equations of the dengue model can be reduced to the following five equations if complete symmetry between the two pathogen strains is assumed: dx/dt = μβ x y 1 -β x y 2 -μ x dy i /dt = β(x + σ z i )y i -(γ + μ)y i (i = 1, 2) dz i /dt = β x y j -σ β z i y i -μ z i (i = 1, 2; j ≠ i).