Confronting the cycle synchronisation paradigm of defoliator outbreaks in space and time—Evidence from two systems in a mixed‐species forest landscape

Defoliators cause extensive damage in boreal and temperate forests of the world. Considerable effort has been invested to understand their individual population dynamics, and despite ample theorising, there is little empirical evidence on factors causing spatial synchrony of pest eruptions at landscape scales. We report on the landscape‐level effect of forest configuration and composition on the intensity of outbreaks of spruce budworm and forest tent caterpillar in a mixedwood boreal forest in northern Minnesota (USA) and adjacent Ontario (Canada), and how this is related to the degree of spatial synchrony in each species' outbreak cycling. Using a large spatiotemporal tree‐ring reconstruction of outbreak impacts across these two systems, we evaluate two contrasting theories governing defoliator outbreaks: harmonic oscillation (a.k.a. ‘clockwork’) and relaxation oscillation (a.k.a. ‘catastrophe’), each with consequences linked to top‐down versus bottom‐up influences on outbreak behaviour in time and space. We find synchrony varies temporally, among outbreak cycles and in direct proportion to cycle peak intensity; however, cycle peak intensities are distributed bimodally in time, and so, therefore, are synchrony coefficients. Spatially, the area where each pest species currently cycles with the greatest peak intensity and synchrony is where their preferred host trees are currently found in greatest proportion. Despite overall synchrony in cycling, we found, in both systems, a persistent negative spatial correlation among successive eruptive pulses of defoliation. Many of these eruptions failed to spread spatially and to coalesce with other spot eruptions to form extensive area‐wide outbreaks. Eruptions often fail to spread at the hardwood‐conifer interface, resulting in outbreak pulses that systematically bounce back and forth between landscape types, particularly when systems were cycling at low amplitude. These over‐dispersed spatial patterns of pulse impact are consistent with a contagious theory of eruption and outbreak spread. They could be considered consistent with harmonic oscillation theory only for populations cycling at different frequencies, with cycling frequency determined by host forest landscape structure. Synthesis. We find that defoliator outbreak dynamics across systems include spatiotemporal signatures of each theoretical paradigm—suggesting a hybrid approach will better characterise outbreak behaviour. Host concentration influences which paradigm dominates the spatial dynamics in any given forest landscape context. Because of the synchronising effect of host concentration on forest insect spatial dynamics, mixedwoods appear to be less prone to intense, synchronised defoliator attack than forests of pure hardwood or pure conifer.


| INTRODUC TI ON
Ecologists have been fascinated with the recurrent periodicity and synchrony of Lepidopteran defoliators for decades (Berryman, 1996;Bjørnstad et al., 1999;Myers, 1988Myers, , 1998;;Myers & Cory, 2013).Yet if forest insect populations cycled both regularly and synchronously, then forest resource managers would have little difficulty predicting the timing and severity of the next forest insect outbreak cycle.Adaptive management to avoid the most severe economic consequences would then be a relatively simple deterministic planning exercise.The harsh reality is that forest pest forecasting remains an elusive objective, particularly under the cascading uncertainties associated with the cumulative effects of forest landscape change and climate change (Dukes et al., 2009).
The growth and durability over the last several decades of the 'clockwork' population cycling paradigm (Bjørnstad & Grenfell, 2001), built upon Christian Huygens' 1673 harmonic oscillation theory (Bell, 1941;Moran, 1953), may seem puzzling if it's contributed so little to our ability to forecast so-called outbreak population 'oscillations'.Part of the reason may be the repugnance of the alternative paradigm-catastrophe theory-which rose to prominence in the 1970s (e.g.Ludwig et al., 1978), prior to the emergence of harmonic oscillation theory under the leadership of Royama (1984Royama ( , 1992Royama ( , 1997)), Berryman (2002) and Turchin (2003).The catastrophe theory of relaxation oscillations suggests that nonlinear stochastic process interactions are the primary driver of forest insect outbreaks, resulting in so-called 'cross-scale dynamics' (sensu Holling, 1973Holling, , 1992;;Raffa et al., 2008).It presents a legitimate concern that well-specified, well-parameterised, specific, causal models of outbreak cycling might be challenging to develop and test if there many key factors that conspire nonlinearly to drive 'cross-scale' outbreak dynamics.
Worse, the critical events that trigger outbreak may occur at population densities far below the limits of detectability using standard operational survey protocols, when populations are endemic, generating only trace indicators of impact.
This paper is about the spatial dynamics of two of the most prominent native periodic forest Lepidoptera in North Americathe spruce budworm (SBW; Choristoneura fumiferana Clem.) and the forest tent caterpillar (FTC; Malacosoma disstria Hbn.)-but in the context of the models used to represent the entire class of forest insect pest to which they belong.SBW and FTC are each univoltine species that time early larval emergence with the spring flush of poorly defended new foliage on their respective preferred hosts, spruce (Picea spp.(Moench) Voss) and fir (Abies balsamea (L.) Mill.) for SBW and trembling aspen (Populus tremuloides Michx.) for FTC.Each exhibit large-scale periodic outbreaks where defoliation of host trees may span millions of square kilometres [(SBW : Berguet et al., 2021;MacLean, 1980) (FTC: Candau et al., 2002;Cooke et al., 2022)].Natural enemies are abundant in both systems and are thought to play an important role in topdown-driven cycling (Roland, 2005;Royama, 1984).However, there is mounting evidence that bottom-up effects from changes in host plant quality or quantity are also implicated in cycling (Nealis, 2016;White, 2018).In the case of SBW, new evidence is emerging on the existence of both positive feedbacks at low-density (i.e.Allee effects associated with mating failure at low density) (Régnière et al., 2013) and host and forest landscape feedbacks at higher densities (Régnière & Nealis, 2019;Robert et al., 2018) that have us revisiting the debate of bottom-up Holling versus topdown Royama (Johns et al., 2019).Indeed, the same is true of the FTC, where we now know that similar low-density Allee effects are present (Evenden et al., 2015), and similar forest landscape feedbacks are involved in regulating outbreak behaviour (Robert et al., 2020).
In our study area, the 'Border Lakes' region of Minnesota and Ontario, Robert et al. (2018) showed that SBW outbreak patterns tend to exhibit greater synchrony in host landscapes dominated by their preferred hosts, while Robert et al. (2020) showed the same for FTC, suggesting that, for both species, bottom-up feedbacks may be a significant factor affecting cycle synchronisation.No attempt was made in these studies to describe the pattern dynamics could be considered consistent with harmonic oscillation theory only for populations cycling at different frequencies, with cycling frequency determined by host forest landscape structure.6. Synthesis.We find that defoliator outbreak dynamics across systems include spatiotemporal signatures of each theoretical paradigm-suggesting a hybrid approach will better characterise outbreak behaviour.Host concentration influences which paradigm dominates the spatial dynamics in any given forest landscape context.
Because of the synchronising effect of host concentration on forest insect spatial dynamics, mixedwoods appear to be less prone to intense, synchronised defoliator attack than forests of pure hardwood or pure conifer.

K E Y W O R D S
Choristoneura fumiferana, cycle amplitude, dendroecology, forest landscape structure, Malacosoma disstria, outbreak dynamics, synchrony of outbreaks at finer spatial and temporal scales.In each case, the focus was on long-term patterns regardless of the causal mechanism of fluctuation.
We revisit those same field data sets and present novel observations that raise important questions about the suitability of the clockwork theory of harmonic oscillation to describe the dynamics of these systems.To understand our line of inquiry, it helps to first review the primary elements of the theory.Appendix S1 was developed to provide a basic summary of the kinds of patterns we expect from simple autoregressive models intended to simulate delayed negative feedback from the top-down action of specialist natural enemies when the cycle-generating process is spatially homogenous in its governing parameters (i.e. a traditional model implementation of the clockwork paradigm; Fleming et al., 1999).
Meanwhile, Box 1 was developed to highlight the contrasting expectations from harmonic oscillation theory versus relaxation oscillation theory.
A key question left unanswered in Robert et al. (2018Robert et al. ( , 2020) ) was the abruptness of change through time in synchronisation dynamics-a hallmark of catastrophe theory that is difficult to reconcile with oscillation theory.
In contrast, according to harmonic oscillation theory, cycles should be of constant amplitude if the stochastic drivers are of constant variance (see Appendix S1), and any changes in cycle amplitude resulting from changes in stochasticity variance should be relatively smooth, not abrupt.Also, cycles of higher amplitude should be better synchronised.According to relaxation oscillation theory, fluctuations must be of discretely varying amplitudes as the system flips abruptly between the lower and upper levels of the outbreak manifold surface.Also, cycles of higher amplitude will only be better synchronised if forest conditions are uniformly supportive of outbreaks.These contrasting expectations are highlighted in Box 1.
In this paper, we seek to determine which of the two outbreak paradigms might be more applicable, and whether the two species systems are behaving similarly.Tree-ring data are recognised as a valuable source of information on the long-term dynamics of outbreak populations causing 50%-100% defoliation but are considered unreliable for making inferences about population dynamics at the low densities around which endemic populations typically transition abruptly to epidemic behaviour, when rates of defoliation are less than 1% (Johns et al., 2019).Consequently, we sought to examine spatiotemporal patterns of outbreak behaviour relative to contrasting predictions from each paradigm, while at the same time avoiding narrowly focused questions on the endemic-epidemic transition that cannot be answered using tree-ring data.Specifically, we attempted to determine: 1. whether cycle synchrony is proportional to cycle amplitude, which would be consistent with a simple linear model of harmonic motion under constant perturbing variance, but inconsistent with a nonlinear model of catastrophic eruption; 2. whether spatial patterns of correlation change smoothly or abruptly in space (e.g. at an ecological boundary), abruptness being more suggestive of a qualitative change in dynamics at some discrete threshold; 3. whether there are abrupt changes in cycle amplitude and synchrony through time, which would be inconsistent with a simple linear model of harmonic motion under constant perturbing variance, but might be more consistent with a discrete model of cuspcatastrophe under bottom-up influence; 4. whether spatial patterns of positive correlation expected at fast time scales switch to spatial patterns of negative correlation at intermediate time scales-a bi-phasic pattern of long-term persistence reported by Cooke and Roland (2023) that is not easily reconciled with harmonic oscillation theory of cycle synchronisation, but could constitute one of Gilmore's (1993) 'catastrophe flags'; and 5. whether there is evidence of asynchronous localised eruption occurring within the broader context of large-scale population cycling, and whether forest context influence the propensity for contagious eruptions to develop into synchronised, area-wide outbreak cycles.
Answers to these questions can suggest which of the outbreak paradigms might be more dominant in these systems, or whether outbreak patterns are a hybrid result of both cyclic and eruptive causal processes, as suggested by Sturtevant et al. (2015).
Ultimately, we seek to explain why forest insect populations generate unpredictably episodic disturbance, despite having a supposedly regular pattern of recurrence, and to determine what role host forest landscape structure might play in governing this dynamic.

BOX 1 Synchronisation theory of forest Lepidopteran outbreaks
We illustrate how different models of cycle causation may result in different modes of synchronisation.To build toward an ultimate contrast between what is expected under two population cycling paradigms, it helps to start simply with model of constant cycling mean and variance, and consider what happens as variance grows, and how synchronisation patterns change depending on whether means and variances transition gradually or switch abruptly (Figure B1).
According to the harmonic oscillation model (Royama, 1992) (top four rows of Figure B1), the amplitude of a pest population cycle, and its likelihood of being spatially synchronised, is proportional to the variance of the perturbing forces affecting the cycle: lowvariance stochasticity generates low-amplitude oscillations (Figure B1a); high-variance stochasticity generates high-amplitude oscillations (Figure B1d) (see Appendix S1 for a demonstration of this phenomenon).Spatial synchronisation emerging from these forces may be driven by random, but spatially autocorrelated, weather effects (Moran, 1953), by dispersal (Barbour, 1990;Royama, 1980), or by the joint action of weather effects on dispersal (Royama, 1984).The degree of spatial covariance in between populations in space may be represented by the relative width of variation around an ensemble of populations depicted in Figure B1.
As a harmonic oscillation, cycle peaks in this model exhibit a normal distribution when the variance of the perturbing force is constant (Figure B1b,e).By extension, we expect a compound distribution (Figure B1h) when the variance of the perturbing force changes gradually over time (Figure B1g).As a result, spatial synchrony may vary through time, but it should do so smoothly, in response to slow changes in environmental variances.As temporal means and variances rise, so would spatial means and variances (contrast Figure B1a,d, and observe the trend in Figure B1g), but spatial covariances would rise even faster, leading to high-amplitude cycle synchronisation (Figure B1i), and a narrower range in observations as variance rises.If the means and variances in a linear oscillation model rise abruptly due to non-stationarity in some model parameters (Figure B1j), then the distribution of cycle peaks will be bimodal (Figure B1k), and the pattern of synchronisation will also be bimodal, rising abruptly with abruptly rising variance (Figure B1l).
Appendix S1 provides a more fulsome and transparent demonstration of these effects in the case of the harmonic oscillation model.
According to the relaxation oscillation, or 'cusp-catastrophe' theory of insect outbreaks (Ludwig et al., 1978), populations may exhibit low-amplitude fluctuations on the lower sheet of an outbreak manifold, but will flip abruptly, synchronously, and catastrophically toward high-amplitude outbreaks when host forest conditions (which generate the 'cusp') permit (Figure B1m).As Gilmore (1993) discusses, in nonlinear systems of this nature, one expects 'catastrophe flags', including: (1) sudden transitions in outbreak behaviour; (2) multi-modal distributions of cycle peaks; and (3) long-term memory in outbreak patterning.Relaxation oscillations produced via nonlinear trophic interactions are not amenable to phase-synchronisation because the cycle phase is not free to be determined by external perturbations, as would be true for a 'phase-forgetting quasi-cycle' (sensu Nisbet & Gurney, 1982), or harmonic oscillation.
Instead, catastrophic outbreaks erupt synchronously only when forest conditions favour high rates of successful reproduction and dispersal, and the outbreak cycle terminates (i.e. the oscillatory impulse 'relaxes') only once the forest is depleted.According to this paradigm, predictability comes only when the critical environmental thresholds favouring outbreak are known with certainty, and when population status is monitored closely for its proximity to the unstable eruption threshold.We expect spatial synchrony under this paradigm only occurs in response to spatially correlated triggering by external factors, as the internal dynamics are so highly nonlinear as to be unpredictable, leading to higher variability and forecasting uncertainty with high-amplitude cycles than with lowamplitude cycles (Figure B1n).
The 'cusp-catastrophe' is the antithesis of the 'clockwork' concept that has inspired so much of the recent literature.Bimodality of cycle amplitudes (Figure B1n) is not necessarily diagnostic of relaxation oscillation, because the same distribution of cycle peaks is observed in an abruptly non-stationary linear oscillation (Figure B1k); however, the abruptness of transitions in cycling amplitude (Figure B1j,m) is not compatible with a low-dimensional linear delayed feedback model of harmonic oscillation, either with constant variance (Figure B1a,d) or drifting variance (Figure B1g).We expect the pattern of synchronisation in the relaxation oscillation (Figure B1o) to be considerably weaker and may even oppose that of the linear oscillation (Figure B1l).This is because external factors, such as the bottom-up effects of host abundance and concentration, are a critical precondition for synchronisation, while dispersal applied to the relaxation oscillation model is most likely to produce travelling waves rather than synchronise disparate populations per se.
The simulations and conceptual schematics of Figure B1 are all based on relatively simple models of cycle causation and synchronisation but serve as a varied template for interpreting observed patterns of outbreaks in spatially replicated long-term observational data.While much of the discussion from theoretical and statistical ecology tends to focus on the simple pattern in the top row (whose origins are demonstrated in Appendix S1), the patterns in the lower rows are more complex, and serve to broaden our perception of the array of possibilities for outbreak species of forest Lepidoptera.aspens (Populus tremuloides Michx., P. grandidentata Michaux)), as well as several species near the northern limit of their range, such as white pine (Pinus strobus L.), red pine (P.resinosa Ait.) and red maple (Acer rubrum L.).Further details and the history of disturbance in the area are described by Robert et al. (2018Robert et al. ( , 2020)).The divergent landscape structures resulting from distinctly different management legacies (Sturtevant et al., 2014) can be readily visualised in a comparative map of satellite-based estimates of forest species composition, which shows a greater concentration of SBW host-tree species at the centre and northeast of the study area, and a greater concentration of FTC host-tree species around the southwestern perimeter of the study area, particularly in Minnesota (Figure 1).

| Study sample design and outbreak reconstruction
We reconstructed outbreak histories of SBW and FTC by sampling tree-rings from 100 mixedwood sites containing 50 white spruce and 50 trembling aspen sample locations, each distributed according to a consistent sampling design throughout the BLL.Field samples were acquired under permits granted by Quetico Provincial Park in Ontario and Lake Superior National Forest in Minnesota.Sample sites were stratified by management legacy (Minnesota managed, Ontario managed, and unmanaged Wilderness in both Minnesota and Ontario), then by longitude (east, centre, west), and latitude (north, south) (Figure 1).These data were aggregated from the level of individual trees to sites to subareas of sites in close proximity (i.e.within 25 km), thereby generating 16 spatially distributed SBW outbreak chronologies, and 15 spatially distributed FTC outbreak chronologies, with a high degree of congruency in the spatial distribution of reconstruction chronologies, despite the overall negative correlation in the coarse-scale distribution of conifers versus hardwoods across the boreal-Laurentian ecotone (Figure 1).This is a landscape for which there are no nonhost tree species that are not defoliated by some other herbivore species and thus might serve as a surrogate 'climatic control' in outbreak reconstruction.However, our reconstruction methods for FTC have been independently validated in Ontario (Cooke & Roland, 2007) and in Minnesota (Cooke et al., 2022), which is a rarity in dendroecology.

F I G U R E B 1
Conceptualised cycling of a model forest tent caterpillar outbreak process, under varying assumptions about cycling mechanics and variances.Log 10 density on the y-axis of the left column of panels represents a relative population density.Black lines in the left column represent mean dynamics of the populations, shaded regions represent conceptualised ranges or confidence intervals for a series ensemble, and open circle symbols represent cycle peaks that are used to define the frequency histogram in the middle column.Blue trend lines are mean expected cycle peaks.Red wiggly lines are the deterministic dynamics underlying the stochastic realisation illustrated (black line).Ellipses within the right column represent hypothetical point clouds illustrating correlations (positive or negative) between cycle peak values and relative synchrony and may be conceptualised as emerging from the variation represented by the shaded regions of the ensemble populations in the left column: the tighter the range on the ensemble, the higher the synchrony.Panels (a) and (d) and histograms (b) and (e) are taken from simulation output from a linear harmonic oscillator (Appendix S1, Figure A2).Histograms based on just 10 cycles tend toward obviously normal distributions, and this tendency becomes increasingly evident as the number of independent cycles (in space or time) increases (see Appendix S1, Figure A1).The remaining panels are conceptual schematics.Rows three and four are logical extensions of rows two and three, given gradually increasing amplitude (g) or an abrupt change in amplitude (j).The bottom row (m) represents expectations of the scenario in row four, but from a nonlinear relaxation oscillator, under which an abrupt transition in amplitude is characteristic.

F I G U R E 1
The spatial distribution of host tree species and sampled stands across the study area.Circles indicate subareas within which tree-level data were aggregated to form outbreak chronologies, 15 for aspen/forest tent caterpillar and 16 for spruce/spruce budworm.Subarea labels follow the naming convention.1st letter: 'M' = Minnesota managed, 'O' = Ontario managed, 'W' = unmanaged wilderness.2nd letter: 'E' = east, 'C' = centre, 'W' = west.3rd letter: 'N' = north, 'S' = south.Groups 1 and 2 were formed on the basis of similarity in forest landscape structure variables (see Section 2 and Figure SI3).
Detailed methods of outbreak reconstruction for SBW and FTC were described by Robert et al. (2018Robert et al. ( , 2020)), respectively.All data transformations and analyses were conducted using the R statistical software package ver.4.2.3 (R Core Team, 2023), as documented in the archived package of code and data (Cooke, 2023).

| Landscape partitioning
The original experimental design in Robert et al. (2012) was based on an assumption of a strong difference in the spatial pattern of forest disturbances across three management zones in Ontario (coarse patchiness), Minnesota (fine patchiness) and the unmanaged Wilderness (limited patchiness) between these two actively managed regions.In subsequent papers, landscape structure variables (both patchiness variables and species composition) were used as independent predictors of local outbreak dynamics of SBW and FTC, respectively (Robert et al., 2018(Robert et al., , 2020)).In this paper, we sought to make the most direct comparison possible between the two defoliator species' responses to contrasting host preferences.We therefore used cluster analysis (function hclust() in R package 'vegan') to group chronologies according to forest landscape variables, ignoring all administrative boundaries.Independent variables consisted of both forest composition and land cover structure variables.In brief, forest composition was derived from Thapa et al. (2020) to represent the relative basal area (reported as a percentage of the total basal area) for both preferred and the full host complement for each respective defoliator, measured at two time periods (1985 and 2005).
Four land cover structure variables (area-weighted mean non-forest patch size (ha), Shannon index of forest cover diversity, percent forest, and forest edge density (m/ha)) were derived from a time-series of land cover maps (Wolter et al., 2012), averaged across the time period 1975-2000.Four landscape structure variables and four species composition variables (measured at two distinct time-points, 1985 and 2005), all measured over neighbourhoods 30 km in radius (yielding a total of 12 forest landscape structure variables), were used to partition the landscape into distinct contiguous groups.The 16 SBW chronologies were derived from 50 distinct sample plots, and the 15 FTC chronologies were derived from another 50 distinct sample plots, yielding 100 plots total for characterising forest landscape structure (Figure 1).
Based on cluster analysis, the BLL clearly could be partitioned into two distinct groups, based on species composition (spruce-fir vs. aspen-hardwood) and structure (% forest vs. edge density) (Table 1).
The forested sites of the Wilderness Area (W) and Ontario East (OE) were more heavily forested, less fragmented and had higher sprucefir content than those of Minnesota (M) and Ontario West (OW) (Figure SI1).MANOVA showed that these differences were highly significant between the two groups (p < 0.001).

| Time-series analysis of aggregated data
To determine whether it was reasonable to use simple autoregressive models to describe forest insect outbreak activity, we conducted a time-series analysis of the aggregated master outbreak a RBA = relative basal area (BA of tree species or grouping divided by the BA; %).
b Hardwood basal areas exclude maple.
c Budworm hosts include white spruce, black spruce, balsam fir.
d Averaged value for the period 1975-2000.

TA B L E 1
Forest landscape attributes for the two groups identified in Figure SI3.

| Distribution of outbreak cycle peaks
We examined the distribution of cycle peaks in the published master outbreak chronologies to determine if they conformed to the simple expectation of a unimodal distribution (see Appendix S1).This was done at the level of the entire study area (as with time-series analysis), but also over the entire set of 16 or 15 subarea chronologies, for the time-interval common to all series (which started in 1928).We were specifically interested in whether de-aggregating the area-wide master chronologies into its constituent subareas altered the distribution of cycle peaks, or whether they maintained a distribution that matched the theoretical expectation of normality set out in Appendix S1.We used a Shapiro-Wilk test to determine if the resulting distribution of cycle peaks conformed to expectations for a normal distribution.
(See Appendix S1 for a demonstration of this expectation.)

| Time-varying synchrony
Cycle peaks having been identified algorithmically, we next defined the troughs between peaks as the years in which outbreak activity was minimum.This in turn allowed us to define an outbreak 'interval' as the time between successive troughs.This was the basis for subsequent analysis of synchrony within and between intervals.A 'cycle' thus refers to the increasing, then decreasing, pattern of pest activity between the troughs defining each interval.This allowed us to distinguish between 'fast' dynamics as patterns of synchrony occurring within an interval, and 'intermediate' dynamics as patterns of synchrony occurring between successive intervals.
For each identified cycle in the master chronology, we computed the mean pairwise correlation between all subarea chronologies ( 16for SBW, 15 for FTC) over all outbreak intervals.We were interested in the possibility that outbreak cycle synchrony was variable across the time-series, and was, in any one interval, related to the amplitude of the oscillation (as measured by the cycle peak intensity) through that interval (our first introduced hypothesis).This would imply a non-stationarity in the global time-series that would seriously violate the assumptions of time-series analysis.Moreover, it would imply that there might be consequences in the spatial domain, in terms of complexity in the spatial patterning of outbreaks, including the possibility of a bi-phasic waxing and waning of patterns of synchronisation.showed that SBW outbreak cycles in the BLL tended to decay to asynchronous more quickly in landscapes that were more hardwood-dominated than conifer-dominated.Robert et al. (2020) used the same method to test whether FTC outbreaks displayed a similar response, yielding the opposite spatial pattern, with more synchronous oscillations occurring in the aspen-dominated patchy landscapes of Ontario and Minnesota, and asynchronous pulses occurring in the less-disturbed, conifer-dominated Wilderness area at the centre of the study area.In this paper, we re-calculated spatial non-parametric covariance functions for just two landscape units that contrasted heavily (see Section 2.3; Figure 1), to provide a direct comparison of relative rates of synchronisation in the two systems.Our expectations were that spatial non-parametric covariance functions for each species would exhibit an abrupt change at the interface between forest landscape groupings (our second introduced hypothesis).

|
Whereas Robert et al. (2018Robert et al. ( , 2020) computed a single spatial non-parametric covariance functions across many decades, one for each landscape zone (and one for each species), we were interested in the possibility that the spatial non-parametric covariance functions changed abruptly through time, as populations alternated between high-amplitude and low-amplitude cycling (our third introduced hypothesis).

| Spatial cross-correlations among outbreak pulses
To measure the abruptness of change in synchronisation dynamics Without knowing the time scale over which an ecological memory process might be operating, it was logical to test for such an effect at additional time scales not rigidly pre-determined by the period of the cycle.We repeated the above analysis, but using a systematically varying number of window frames, or intervals, from 3 to 13, and with the interval boundaries determined randomly, through a bootstrap resampling approach.We expected that spatial correlations in spatial impacts between successive intervals would become increasingly negative the slower the time scale of analysis, consistent with the same long-term persistence pattern identified by Cooke and Roland (2023).The study-area wide outbreak chronologies were short relative to the cycling frequency of the outbreak processes, yielding just 6 and 10 cycles respectively for SBW and FTC (Figure 2a,g).The cycle peak intensities were highly variable, leading to non-normal distributions; however, they were also characteristically noisy (Figure 2b,h).
Consequently, we further characterised the distribution of SBW and FTC cycle peak intensities measured at the finer level of the 16 or 15 subarea chronologies (Figures SI1 and SI2).We found cycle peak distributions at this scale were not normally distributed (Shapiro-Wilk p < 0.0001), but rather bimodal, with SBW cycles that peaked at intensities of 35% and 95% of trees affected (Figure 3a), and FTC cycles that peaked at intensities of 15% and 65% of trees affected (Figure 3b).
We also examined cycle peak distributions between forest landscape groupings.Although cycle peak intensities for each herbivore species appeared more bimodal where their respective host was less concentrated, a Kolmogorov-Smirnoff test indicated differences were not significant (Figure SI2).The distributions of cycle peak intensities were generally bimodal at all scales and over both forest landscape groupings.

| Time-varying synchrony
For both defoliator species, outbreak synchrony, r i , varied significantly among outbreak intervals (Figure 2a,g), and was positively correlated with the cycle peak intensity during each interval.For SBW, the correlation between interval synchrony among 16 chronologies and mean outbreak cycle peak intensity was 0.99 across the five identified outbreak cycles for which a correlation could be measured (Figure 2f).For FTC, the correlation between interval synchrony among 15 chronologies and mean outbreak cycle peak intensity was 0.82 across the 10 identified outbreak cycles (Figure 2l).For both defoliator species, more than half the outbreak cycles identified peaked at an intensity lower than 50%.For SBW, outbreak cycles II-IV were all low-intensity events where the interval synchrony did not exceed 0.26, constituting 91 years of low synchrony (i.e.1887-1978, Figure 2a), which is more than half the length of the 1836-2006 chronology.For FTC, outbreak cycles II-IV and VII-VIII were all low-intensity events where the interval synchrony did not exceed 0.25, constituting 66 years of low synchrony (i.e.1888-1931and 1961-1984, Figure , Figure 2g), which is half the length of the 1872-2006 chronology.In both systems, the total amount of time spent cycling at low amplitude with low synchrony was roughly 50%, with stretches that endured roughly two to three cycles.In both systems, the last outbreak cycle (cycle VI in the 1980s for SBW, and cycle X in the early 2000s for FTC) represented events with the highest intensity and synchrony for each respective tree-ring record.

| Host-dependent patterns of synchrony decay
When we computed SNCFs for SBW and FTC outbreak chronologies according to the empirically defined landscape groupings (Figure 1; Figure SI4), we found the two insect species exhibited opposing patterns of synchronisation, with SBW outbreak cycles more synchronised in the heavily forested, conifer-dominated Group 2 (Figure 4) and FTC outbreak cycles more synchronised in heavily disturbed, aspen-dominated Group 1 (Figure 5).A greater quantity of a given host group tends to favour a higher level of outbreak cycle synchronisation for herbivores exploiting that host group.
When we recomputed spatial non-parametric covariance functions to accommodate variable synchrony among outbreak intervals, the landscape response became even clearer.For SBW, the interval synchrony rose from cycle IV to V to VI, in both landscape groupings 1 and 2, and synchrony was always higher within the conifer dominated group 2 plots (Figure 4).SBW outbreak covariances dropped below zero in the aspen-dominated group 1 plots during cycle IV.For FTC, the interval synchrony was lowest during The distribution of cycle peaks for spruce budworm (SBW) (a) and forest tent caterpillar (FTC) (b) when ampd() is applied to the 16 SBW and 15 FTC subarea chronologies (which are illustrated in Figures SI1 and SI2).The distributions are clearly bimodal, with the distribution mean represented as a vertical line.Shapiro-Wilk test for normality is rejected at p < 0.001 in both cases.Although the mean cycle intensity is roughly 50% for both species, the mean is not a good descriptor of central tendency in either case due to the bimodality.SBW outbreaks at the subarea level tend to be of moderate or severe intensity (35% vs. 95% of trees in a chronology affected).FTC outbreaks tend to be of light or moderate intensity (15% vs. 65% of trees in a chronology affected).These distributions may be compared to those in Figure 2b,h for study-area wide chronologies, or contrasted against those in Figure A1 for theoretical expectations from a secondorder autoregressive cycling model.
the middle cycles VII, VIII, and was higher during the earliest and latest cycles, VI and X (Figure 5).During the low-amplitude FTC cycles VII-IX, synchrony did not differ between landscape groupings 1 and 2. During the high-amplitude cycles VI and X, the pattern switched, with synchrony being higher in group 2 during VI and in group 1 during X.As with SBW, FTC outbreak covariances sometimes dropped below zero during the low-amplitude cycles VII to IX-whenever the system was cycling with low amplitude.

| Spatial asynchrony
When we plotted the spatial pattern of growth impacts among SBW outbreak cycles we observed negative correlations between successive intervals (Figure 6).Calculating spatial correlations among the two pairs of intervals, we found that the correlation between successive intervals averaged r = −0.23 for the 16 subareas and −0.15 for the 50 sites.When we plotted the spatial pattern of growth impacts among FTC outbreak cycles, the correlation between successive intervals was similar, averaging r = −0.23 for the 15 subareas and −0.11 for the 50 sites (Figure 7).Each herbivore showed a similar pattern: outbreak activity was concentrated within different places through time in a somewhat systematic fashion.
When we relaxed the constraint on the number and position of windowed intervals and explored a wide range of windows, we found that the spatial correlations in impacts between successive intervals dropped considerably, and the slower the time scale of analysis, the stronger the repulsion between successive outbreak pulses (Figure SI3).

| Host effect on contagious partial cycling
For both species, fluctuations at the fastest time scale of the subinterval and the finest spatial scale of the sub-area were complex, F I G U R E 4 Spatial non-parametric covariance functions for spruce budworm during cycle intervals IV, V, and VI, and the entire period IV-VI for aspen-dominated Group 1 (M+OW; black line, heavy-shaded 95% CI) versus conifer-dominated group 2 (W+OE; red, light-shaded 95% CI). Figure SI3 maps these groupings.
and not consistent with a simple hypothesis of spatially synchronised cycling.While algorithm ampd() identified three and six cycle peaks for SBW and FTC, respectively, at the level of the whole study area (Figures 6 and 7), additional local peaks could be identified at the subarea-level (see Figures SI4 for SBW and SI5 for FTC).VMRs for impact maps at this scale were heavily skewed to values much larger than one, indicating highly over-dispersed, not diffuse, patterns of impact (Figure 8).In fact, VMRs on site-level chronologies were biased so high that the SNCFs behaved erratically at short distance classes.[For SBW and FTC examples,see Figures SI6 and SI7,respectively.]The confidence intervals on correlograms were extremely wide at short distance classes, particularly for SBW, including both positive and negative values-indicating that heavy impacts at times were restricted to a given site, while low impacts were occurring at adjacent sites.The spatial pattern of impact thus frequently included both contagious eruption and local repulsion.This furthermore explains how spatial correlations between intervals rose from <−0.2 to >−0.2 as the resolution of analysis was raised from 15 or 16 subareas to 50 sites.
Among pulses, spatial patterns of impact would sometimes correlate positively with forest attribute principal components, and sometimes negatively (Figure 9).In the case of SBW, the distribution of forest-impact correlations was bimodal (Figure 9b,c), indicating a tension in eruptive cycling between the two forest landscape units-a tension that was strong enough to resist the many forces of cycle synchronisation that should be operating at that scale, from dispersal effects to weather effects (see Appendix S1).For FTC, the distribution of forest-impact correlations was just as wide-including both positive and negative values-but also including several neutral values for forest principal component 1 (Figure 9d).The motion in pulses among forest landscape units was thus not as abrupt for FTC as for SBW.For both systems we observed cycles that were somewhat synchronous across large parts of the study area, but these were frequently associated with localised pulses of impact that tended to shift

| Irregular cycling
Although defoliation impacts in both the SBW and FTC systems tended to be periodic over time, with major outbreaks of SBW occurring every 45 years and major outbreaks of FTC occurring every 15 years, the oscillations in each were irregular.Cycle amplitude in both systems tended to vary slowly, with extended periods of low amplitude cycling enduring over decades (Figure 2).During this time, spatial synchrony was low, and the resulting spatial patterns of impact indicated outbreaks that consistently failed to expand to impact every plot (Figures 6 and 7).In both systems, the contrast between high-versus low-amplitude cycling was sharp, with cycle peak distributions that were bimodal at all spatial scales, ranging from the subarea (~20 km) up to the entire study area (~200 km) (Figure 3;

| Variable synchrony
The observation that synchrony in both systems was higher when cycling was high-amplitude is consistent with the modelling premises of Barbour (1990) under the harmonic oscillation paradigm (Box 1).However, a sharp discontinuity between high-amplitude and low-amplitude cycling appears to be a new observationmore consistent with a relaxation oscillation model (Box 1).(1984) argued that density-dependent dispersal tends to raise cycle amplitudes in simple autoregressive feedback models, and Régnière and Nealis (2019) demonstrated SBW dispersal is indeed density-dependent.Perhaps the bimodal distribution of cycle peaks is a simple manifestation of a positive feedback surrounding density-dependent dispersal: when a system is cycling with lower amplitude, it produces fewer dispersers, further hampering the system's propensity for synchronisation.Other, more F I G U R E 6 Spatial patterns of the impact of spruce budworm on white spruce growth across the study area during three windows of spruce budworm outbreak cycling.The average spatial crosscorrelation between successive outbreak frames is −0.23 for 16 subareas and −0.15 for 50 sites, indicating a systematic motion in spatial patterning of peak impacts.Spatial interpolation among sample locations was done using inverse distance weighted regression surface (weighting power = 0.35).The expected correlation is 0.00 (Appendix S1, Figure A7a).

Royama
elaborate mechanisms could be involved, such as the operation of a slower ecosystem process arising in either the natural enemy community or through landscape-level host-plant effects.These possibilities will be discussed.
In both systems, defoliation impacts in any given year tended to be positively spatially autocorrelated; however, these positive spatial autocorrelations tended to decline in space, and, consistent with our first hypothesis, they decline more quickly in landscapes that had a lower proportion of relevant host-tree species, spruce-fir and continuous forest for SBW (Figure 4), aspen and high-edge forests for FTC (Figure 5).

| Comparison to related studies
Our analyses showed that outbreaks during any one cycle interval will often fail to spread to cover a significant part of the landscape, with many cycles affecting just 15% or 35% of trees annually at their peak.Two major SBW cycles were identified through the interval 1830-2005, along with three intermediate cycles and one low-intensity cycle (Figure 2a); however, in many individual series we also observed additional peaks of lesser amplitude occurring between the major peaks (Figure SI4).This is a pattern that has been reported in the SBW dendroecology literature, but not extensively discussed.Bouchard et al. (2006) remarked on the occurrence of SBW cycles in Quebec that failed to impact a majority of spruce and fir in the landscape.Different authors might not consider these low intensity impacts as representing full 'outbreaks'.Fraver et al. (2007), for example, studying white spruce in a small area in Maine, did not count such low-intensity fluctuations as 'outbreaks', and so computed an outbreak cycle frequency of 67 years for SBW, which is half the frequency reported for nearby New Brunswick and Quebec (Jardon et al., 2003;Royama, 1992).
This remarkably low frequency is contingent on the definition of F I G U R E 7 Spatial patterns of the impact of forest tent caterpillar (FTC) on trembling aspen growth reductions across the study area during six windows of FTC outbreak cycling.The average spatial cross-correlation between successive outbreak frames is −0.23 for 15 subareas and −0.11 for 50 sites, indicating a systematic motion in spatial patterning of peak impacts.Spatial interpolation as in Figure 6.The expected correlation is 0.00 (Appendix S1, Figure A7a).

F I G U R E 8
Variance-to-mean ratio (VMR) distributions for sub-interval impact maps, for (a) spruce budworm (SBW) and (b) forest tent caterpillar (FTC).VMR ≫ 1 indicates a highly clustered, over-dispersed pattern of spatial variation in impact.VMR ~ 1 indicates a random spatial distribution of impact.(Spatial patterns of subinterval impact are illustrated, along with VMR computations, in Figures SI6 and SI7 for SBW and FTC.)These distributions may be compared directly to the theoretical expectations from a harmonic oscillation model, where the distribution mean is roughly 1.6 (Appendix S1, Figure A7b).what constitutes an outbreak-inclusion of low-intensity cycles restores the outbreak cycle frequency to 34 years.In that regard, our results concur with those of Fraver et al. (2007), in that we observed a wide range in outbreak cycle intensities.The distribution of outbreak cycle peaks at the scale of the BLL appears to vary between low and moderate intensities (Figure 3).A similar pattern of frequent minor peaks was observed for FTC (Figure SI5).This is consistent with the observation in the literature that mapped FTC defoliation during an outbreak interval rarely covers more than 50% of the host landscape available for defoliation (Cooke et al., 2009).
Our observation that FTC outbreak cycle intensity and synchrony is correlated with host forest cover would appear to clash with other conclusions from long-term studies of FTC population numbers and defoliation patterns in landscapes heavily fragmented by agriculture and land clearing, which show that rather than resulting in shorter outbreaks, forest fragmentation results in more frequent cycling (Roland, 2005), with outbreaks that can last several years longer than in less disturbed landscapes (Roland, 1993).
This contrast highlights the importance of clarifying the landscape context of any study.Our study in the BLL is situated between the Fort Frances and Thunder Bay districts studied by Roland (1993), but far away from active land development associated with those regional centres.It is possible that the modest forest fragmentation in our study, because it is largely a result of lakes and nonhost species forest, is sufficient to affect the level of defoliator inter-patch dispersal, and possibly the species composition of the natural enemy community, but is insufficient to completely impede the movement of defoliator-hunting parasitoids that are critical for precipitating the downward phase of the population cycle (Roland & Taylor, 1997).Fragmented mixedwoods may be more prone to periods of low-intensity outbreak behaviour that manifest as asynchronous partial outbreaks, relative to the classic oscillations that are a standard focus in the theoretical spatial population dynamics literature.
We suggest that bimodally distributed cycle peak intensities may be a relatively common feature in tree-ring outbreak data that has been under-reported, as we see similar patterns of bimodality of impacts of western SBW in British Columbia, for example (Axelson et al., 2015).We encourage the community of landscape dendroecologists to examine their data at multiple spatial and temporal scales to determine if they are recording similar spatiotemporal complexity that could shed additional light on the clockwork-catastrophe debate.
We are not the first to report variable rates of synchronisation over time or through space (Bjørnstad et al., 2008).The observation has been discussed extensively in the case of the invasive Lymantria dispar dispar (spongy moth) in the United States (Johnson et al., 2006) and around the world (Johnson et al., 2005), and the indigenous larch budmoth (Zeiraphera diniana) in the European Alps (Bjørnstad et al., 2002;Büntgen et al., 2009;Hartl-Meier et al., 2017;Johnson et al., 2004;Konter et al., 2015;Saulnier et al., 2017).what can be expected in a highly connected, host-rich landscape.

| Driving mechanism
If forest insect population cycles 'yearn to align' (Bjørnstad, 2000), then what inhibits them from synchronising to the point that all populations peak simultaneously everywhere?Our analysis suggests that when the FTC system is cycling with low-amplitude, then populations in the aspen-dominated group 1 cycle completely out of phase with populations in the conifer-dominated group 2 (Figure SI8).This result is consistent with Barbour's (1990) contention that cycle amplitude is a key ingredient for synchronisation, higher amplitudes resulting in greater numbers of betweenpatch migrants, and thus greater levels of cycle synchronisation.
However, this alone cannot explain the patterns observed, for why are there negative spatial correlations between outbreak pulses that persist for decades?
The simplest possible explanation is that insects of a given species are attempting to cycle at different frequencies in the different landscape types, and this creates tension in the network of harmonic oscillations.Our neutral model (Appendix S1) was premised on a homogenous landscape in which all populations cycled at the same frequency.This assumption may not realistic.If, as reported by Roland (2005), populations tend to cycle at a lower amplitude and faster frequency in host-poor landscapes, this could possibly generate the observed asynchrony in cycling between forest landscape units.Cooke and Roland (2023) reported exactly this kind of harmonic tension at the forest-farmland interface for FTC in central Alberta.
Another possibility is that there is something yet unidentified that systematically suppresses the growth of population cycles, so that outbreaks tend to shift about spatially from one cycle to the next.This is particularly puzzling for a species such as FTC, which is known to cause tree mortality in only exceptional circumstances (Schowalter, 2017).One possibility is that there are other agents operating in the system whose roles are neglected in field population studies.This could include soil-borne pathogens capable of epizootics, or higher-level trophic interactions, such as from hyper-parasitism or even hyper-hyper-parasitism, which does occur under the broad umbrella of SBW community ecology (Eveleigh et al., 2007).Burden et al. (2003) suggested covert infections as a mechanism for long-term persistence of baculoviruses.Meanwhile, Cooper et al. (2003) hypothesised that latent infections of nucleopolyhedrovirus (NPV) may play a role in the persistence dynamics of outbreaks of western and FTCs, including the potential for cross-infectivity.
Notably, both SBW and FTC are subject to epizootic infection by microsporidian species of Nosema; however, their role in shaping cycle intensity is poorly understood.Any of these candidates could serve as the source of ecological memory that counteracts the cycle synchronisation process and limits the extent of outbreaks.

| Confronting the paradigms
Relaxation oscillation theory was inspired by Réné Thom's (1972) catastrophe theory, emphasising the impossibility of precisely predicting the timing of crisis events.Harmonic oscillation theory-advanced by Royama (1984) in the case of SBW and generalised to a broader array of animal systems by Royama (1992)-was inspired by Lotka-Volterra, emphasising the regularity of recurring population cycling events (Wangersky, 1978).The two theories have been contrasted at a high level in terms of their ability to explain population cycling in SBW (Pureswaran et al., 2016), but considerably less attention has been given to how they might relate to observed patterns of synchronisation in actual systems (see Box 1).
The SBW outbreak cycle has been represented mathematically as both a nonlinear relaxation oscillation (Jones, 1974;Ludwig et al., 1978), and as a linear harmonic oscillation (Fleming et al., 1999), as has the FTC, for example, Rose and Harmsen (1981) versus Cobbold et al. (2005) and Roland (2005).Support in favour of each theory has waxed and waned over the decades, with the question remaining unresolved due to the tremendous amounts of data required to disprove either one (Sturtevant et al., 2015).
Harmonic oscillation theory has yet to be revisited in the rest of forest Lepidopteran population ecology, which continues to favour the popular view of the last 20 years-of forest insect outbreaks as a form of 'noisy clockwork' (Bjørnstad & Grenfell, 2001;Haynes & Walter, 2022).However, our investigation suggests it may be time for a third hypothesis: that these clockwork forest pest catastrophes might be regulated jointly by both top-down and bottom-up effects, with the balance of evidence in any study being context-dependent (Kneeshaw et al., 2015;Sturtevant et al., 2015Sturtevant et al., , 2023)).These systems each exhibit some of the 'catastrophe flags' of Gilmore (1993), yet the periodicity of outbreaks occurring at all spatial scales is undeniable.

| CON CLUS IONS
We demonstrated that the SBW cycles most synchronously in areas where its spruce-fir host was most abundant.We simultaneously demonstrated that an unrelated herbivore, the FTC, in the same mixedwood landscape, exhibits an identical response to the distribution of its broadleaf host trees, leading to a compelling result where the two insects are cycling synchronously in opposing areas of the same landscape.
We showed that neither SBW nor FTC cycles sinusoidally, and neither fluctuates in a particularly well-synchronised manner.
Outbreaks frequently fail to rise to an intensity, or spread to an extent, where they are affecting more than half of the host trees available for defoliation.This creates a challenge for proponents of both the cycling and eruptive theories of periodic forest insect outbreak.If one imagines periodic outbreaks are the result of a top-down-driven harmonic oscillation, then why do some oscillations fail to rise to an intensity that affects more than a small fraction of host trees?Conversely, if one imagines forest pest catastrophes are the result of a bottom-up-driven relaxation oscillation, then why do so many contagious eruptions fail to spread to cover the available range of host forest?We found evidence in favour of each hypothesis with consistencies in patterns across both systems.A logical conclusion is that the governing mechanism is hybrid cyclic-eruptive.pattern (Bjørnstad & Grenfell, 2001), we note that it tends to be the aberrations, rather than the regular patterns, that create the most significant forest declines in need of prediction (e.g.Cooke et al., 2022).
For nearly a century, ecologists have debated the role of the forest in shaping insect outbreak dynamics (Kneeshaw et al., 2021).In the case of defoliators characterised by their apparent periodicity and synchrony, there has been considerable scepticism that the forest condition will have much influence on the dynamics of defoliators beyond the simplest case that outbreaks do not occur where there is no host to support them (Miller & Rusnock, 1993;Muzika & Liebhold, 2000).We showed that host forest cover is correlated with both the intensity and the spatial scale of synchrony of periodic outbreaks consistently across two contrasting insect species.This result supports the notion that forest management can play an important role in the resistance and resilience of managed forests affected by insect disturbance, and likely applies to other forest-pest systems, from invasion ecology to indigenous pests from around the world.
Stig Larsson (1989) wrote that it was 'stressful times for the plant stress hypothesis'.We are similarly compelled to ask in 2023 whether 'the clock is ticking' on Huygens' 1673 Horologium Oscillatorium (Bell, 1941)-the ultimate source of the idea that periodic insect outbreaks are the result of the clockwork action of natural enemies combined with a 'Moran effect' (Moran, 1953) that serves to synchronise independently oscillating pest populations.
Our research suggests that the only path available to rescue harmonic oscillation theory from the jaws of catastrophe theory is to chronologies presented inRobert et al. (2018) (1832-2005 for   SBW)  andRobert et al. (2020) (1878-2006 for FTC).We computed correlograms to diagnose autocorrelations and partial autocorrelations and performed spectral analysis to diagnose the level and nature of any periodic behaviour.The functions acf(), pacf() and spectrum() are members of the 'forecast' package in R(Hyndman & Khandakar, 2008).

Function
ampd() (Automatic Peak Detection in Noisy Periodic and Quasi-Periodic Signals) in R package 'AMPD' ver.0.2.(Scholkmann et al., 2012) was used to infer cycle peaks in the field data.The function takes a single scaling parameter, L, to define its sensitivity.For all SBW outbreak chronologies, L was set to L = 8.For FTC outbreak data, where the system cycles much more quickly, it was set to L = 4.The fine tuning of this parameter is somewhat arbitrary, however use of ampd() is preferable to traditional methods in dendroecology, as it removes subjectivity from the definition of an outbreak threshold, thus enforcing a consistent application of the filter throughout a time-series and allowing repeatability in applications across scientific studies.
though time we sought to determine whether spatial patterns of positive correlation expected at fast time scales would give way to spatial patterns of negative correlation at intermediate time scales (our fourth introduced hypothesis).For each of the cycle peaks identified in Section 2.3, we defined the outbreak cycle interval as starting and stopping when the outbreak intensity attained its minimum value between cycle peaks.We then mapped the average percentage of host trees impacted by each herbivore species during each window interval.Spatial correlations were then computed between all pairs of windows, to determine if successive outbreak pulses in each species were occurring in the same area cycle after cycle or if they were moving around in a pattern reflective of the previous outbreak, indicating some form of ecological memory.

F
Time-series analysis of the full aggregate spruce budworm (SBW) (a) and forest tent caterpillar (FTC) (g) outbreak chronologies (thick lines), with cycle peaks inferred using ampd(), histograms of which are shown in (b) and (h).Also plotted in (a) and (g) are the individual subarea chronologies (thin grey lines).These are the basis for the estimates of interval synchrony r i , (square labels at mid-interval) measured as the mean pairwise correlation between chronologies during each outbreak interval (intervals depicted as vertical black bars).Autocorrelation functions in (c), (i) indicate half-cycles of 16 years and 5 years for SBW and FTC, respectively.Partial autocorrelation functions in (d), (j) indicate highly significant positive first order and negative second order feedbacks for both SBW and FTC.Spectral analyses in (e), (k) indicate strong periodicity of 45 and 15 years per cycle for SBW and FTC respectively (series smoothed using Daniell smoother with spans = 5).The among-interval correlation between cycle peak intensities (open circles in (a) and (g)) and interval synchrony, r i , is strongly positive for both defoliator species ((f) and (l); linear model fit shown as red line).Cycle peaks in (a) and (g) are plotted as open or filled black circles depending on whether the number of source chronologies is sufficient to estimate a spatial covariance function.Early intervals (black circles) lack sample depth.This is also reflected in the shading scheme in (f) and (l).Between-interval spatial correlations are indicated in circles at interval boundaries.Mean time-series correlation within intervals is 0.35 for SBW and 0.30 for FTC.Mean spatial cross correlation between intervals is −0.15 for SBW and −0.11 for FTC.
Figure B1), but are somewhat more reminiscent of nonlinear, bottom-up-driven relaxation oscillator theory (bottom row of Box 1, Figure B1).In fact, a careful review of each result in turn suggests a striking similarity to the hybrid patterns illustrated in the fourth row of Box 1, Figure B1.

Figure
Figure SI2).A bimodal distribution of cycle peaks is inconsistent with the unimodal, indeed, normal, pattern expected from the sort of simple autoregressive feedback models that have been used in numerous publications to represent the population dynamics these defoliator species (see Appendix S1).In other words, neither system conforms to a simple theory of harmonic oscillation.Such abrupt changes in cycle amplitude, however, are consistent with the nonlinear relaxation oscillation model (see Box 1).

F
I G U R E 9 (a) Principal components biplot of the 12 forest species composition and structure variables used to define landscape units in Figure 1 and Figure SI4.Note the strong separation of the 100 plots into groups based on species relative basal area (SF = spruce-fir, F = fir, A = aspen, H-M = hardwood less maple), and patch structure (PF = percent forest, AWMPS = area-weighted mean patch size.ED = forestnon-forest edge density, CD=Shannon's index of cover diversity).Species composition was estimated in both 1985 (85) and 2005 (05).All estimates are for neighbourhoods 30 km in diameter.Symbols nudged slightly for enhanced legibility.Histograms in (b)-(e) describe the distribution of correlations between subinterval-scale impact maps and the PC1 and PC2 forest attributes plotted in (a).(b) and (c) are for spruce budworm on PC1 and PC2.(d) and (e) are for forest tent caterpillar on PC1 and PC2.Dashed blue vertical lines indicate splits between modes of the distribution of forest-impact correlations.Multi-modal distributions are evidence that outbreak cycle synchrony tends to break down along forest unit boundaries, as eruptive pulses bounce back and forth between forest landscape units.
each case, weather, mesoclimate, macroclimate and/or climate change have been implicated as potential explanations for the occasional gain or loss of periodicity and synchrony through time.The idea that forest landscape structure, as well as weather and climate change, also might influence spatial synchrony and changes in synchrony through time has been discussed for several species, and arguments have been made that reduced host forest landscape connectivity is associated with reduced periodicity, reduced cycle amplitude and/or reduced cycle synchrony.However, earlier studies suffer from some drawbacks.For larch budmoth, Johnson et al. (2004) did not use independently derived forest structure data to support their contention that fragmentation of the larch forest was driving a breakdown in synchrony of travelling waves.For spongy moth, Haynes et al. (2013) concluded that weather effects on masting, not forest composition per se, drove spatial synchrony in outbreaks.Haynes et al. (2019) reported, contrarily, that forests with abundant oak and hickory host trees were more likely to exhibit the high-frequency 4-5-year cycling that Johnson et al. (2006) associated with pulsed invasion; however, this effect of forest composition on subharmonic cycling is independent of the question of what drives synchronisation in the primary 8-10-year outbreak cycle.In our systems, which are not alien invasive, and therefore do not suffer from the complication of having a moving invasion front that pulses forward as populations surge upward every 4-5 years, we saw no evidence of well-synchronised, high-frequency, subharmonic oscillations being associated with the primary host forest cover type.The asynchrony we observed was between areas where populations were cycling at similar frequencies to the dominant frequency: 45 years for SBW, and 15 years for FTC.With SBW, it was the cycle amplitudes that were variable during the era of low cycle amplitude and low synchrony (FigureSI9).With FTC, it was the cycle phases that failed to align during the era of low cycle amplitude and low synchrony (FigureSI8).In each case, partial autocorrelation coefficients that were weakly negative at lags of 3-6 years at the subarea level (averaged across subareas in FigureSI10) accumulated over those short time scales to induce a long-term persistence effect, or memory effect, whereby the occurrence of a cycle in a given area tended to attenuate the next cycle in that area.Our results nevertheless support the broad interpretations ofLiebhold et al. (2022) in spongy moth andJohnson et al. (2004) in larch budmoth: when systems lose periodicity and amplitude, they lose synchrony, and the loss of synchrony is greater in heterogeneous landscapes than in host-rich landscapes.The reason our results provide valuable support for their work is because we are studying forest insect outbreak dynamics in the relatively intact forest of the Border Lakes region of North America, which provides us a solid baseline for quantifying Regardless which theory one favours, there is a demonstrable role of host forest abundance and distribution in regulating the spatial dynamics of synchronisation of eruptive cycles.Increasing host forest cover tends to attenuate asynchronous eruptive behaviour, promoting synchronous cycling.The fact that SBW exploits sprucefir and FTC aspen means that the two insects tend to cycle more synchronously in opposing parts of the study landscape-the FTC in the aspen-rich, edge-dominated forests of Minnesota, and the SBW in the more heavily forested, conifer-rich forests of the Border Lakes Wilderness region.We found forested areas not impacted during a given outbreak cycle interval were the ones most likely to be impacted in the next interval.Such ecological memory is a novel result that contradicts both the most frequently cited theories of forest insect population dynamics (Appendix S1) and common practice in operational pest hazard analysis in forestry, where positive spatial and temporal autocorrelations are tacitly assumed: what one observed spatially during an earlier window in time is what one is most likely observe during the next window.The opposite appears true, at least at the scale of this study.While one might attribute such complexity to 'noise' around a relatively robust clockwork-like outbreak suppose that (a) different populations in different landscapes cycle at different frequencies; (b) the factors that modulate cycle amplitude can switch modes for decades at a time, such that sustained low-amplitude cycling is commonplace, resulting in (c) insufficient migration to reconcile spatial differences in cycling frequency.Even then, there is a clear role for bottom-up effects of host forest in synchronising cycles, so the cycle-generating process cannot possibly be second-order autoregressive.Even if faster 'top-down' forces (i.e.predator-prey-type interactions) are responsible for the induction of population cycling, these cycles are subject to additional regulation by slower 'bottom-up' effects of host-plants exerted at the landscape level, as outlined bySturtevant et al. (2015).