Forest landscape structure in ﬂ uences the cyclic-eruptive spatial dynamics of forest tent caterpillar outbreaks

. A fundamental question in forest insect ecology is the role of forest landscape structure, partic-ularly the amount and spatial con ﬁ guration of host tree species, in shaping the dynamics of recurring forest insect outbreaks. For forest tent caterpillar (FTC), independent studies do not converge on a singular conclusion, although all indicate that forest structure in ﬂ uences outbreak dynamics. These studies also vary in how they treat climate as a covariate. We evaluated the relative importance of host forest landscape composition and con ﬁ guration, as well as climate, for their in ﬂ uence on FTC outbreak cycling in the twen-tieth century. We predicted that FTC outbreaks would exhibit greater synchrony and intensity within areas associated with higher abundance of host trees. We reconstructed FTC outbreaks from 1928 to 2006 using tree-ring analysis within a well-structured experimental landscape located in northwestern Ontario and northern Minnesota. Time-series clustering and spatial nonparametric covariance were used to determine whether similarities in time series and patterns in spatial synchrony corresponded with land management history. Using constrained ordination, we compared statistical properties of outbreak time series to landscape variables representing host abundance, forest con ﬁ guration, and climate. We found no evidence of climatic effects at the scale of this study, but a signi ﬁ cant albeit small in ﬂ uence of landscape structure on outbreak dynamics. Outbreaks were more synchronous and more cyclic within managed zones containing a greater relative abundance of aspen and other hardwood host tree species, compared with the more coni-fer-dominated Wilderness area. Yet, we also observed asynchronous outbreak dynamics across the study area, such that correlations with slower-changing forest landscape variables varied starkly among outbreak pulses. Consequently, the strength of relationship between landscape variables and FTC outbreak patterns varied strongly through time — a result that may explain why short-term studies yield conclusions that are at odds with one another. Our results speak to the importance of long time series, contrasting landscape structure, use of multivariate methods, and controlling for climatic variation when investigating the effects of forest landscape structure on the cyclic-eruptive spatial dynamics for forest defoliators.


Abstract.
A fundamental question in forest insect ecology is the role of forest landscape structure, particularly the amount and spatial configuration of host tree species, in shaping the dynamics of recurring forest insect outbreaks. For forest tent caterpillar (FTC), independent studies do not converge on a singular conclusion, although all indicate that forest structure influences outbreak dynamics. These studies also vary in how they treat climate as a covariate. We evaluated the relative importance of host forest landscape composition and configuration, as well as climate, for their influence on FTC outbreak cycling in the twentieth century. We predicted that FTC outbreaks would exhibit greater synchrony and intensity within areas associated with higher abundance of host trees. We reconstructed FTC outbreaks from 1928 to 2006 using tree-ring analysis within a well-structured experimental landscape located in northwestern Ontario and northern Minnesota. Time-series clustering and spatial nonparametric covariance were used to determine whether similarities in time series and patterns in spatial synchrony corresponded with land management history. Using constrained ordination, we compared statistical properties of outbreak time series to landscape variables representing host abundance, forest configuration, and climate. We found no evidence of climatic effects at the scale of this study, but a significant albeit small influence of landscape structure on outbreak dynamics. Outbreaks were more synchronous and more cyclic within managed zones containing a greater relative abundance of aspen and other hardwood host tree species, compared with the more conifer-dominated Wilderness area. Yet, we also observed asynchronous outbreak dynamics across the study area, such that correlations with slower-changing forest landscape variables varied starkly among outbreak pulses. Consequently, the strength of relationship between landscape variables and FTC outbreak patterns varied strongly through time-a result that may explain why short-term studies yield conclusions that are at odds with one another. Our results speak to the importance of long time series, contrasting landscape structure, use of multivariate methods, and controlling for climatic variation when investigating the effects of forest landscape structure on the cyclic-eruptive spatial dynamics for forest defoliators.

INTRODUCTION
The forest tent caterpillar (FTC; Malacosoma disstria H€ ubner) is one of the most widely distributed forest insects in North America (Stehr and Cook 1968), exhibiting periodic outbreaks that span millions of square kilometers in eastern Canada (Cooke et al. 2012), and impacting forests as far south as the Gulf Coast of the southern United States (Parry and Goyer 2004). These outbreaks occasionally cause the decline of millions of hectares of aspen (Populus spp.) (Churchill et al. 1964, Candau et al. 2002 and sugar maple (Acer saccharum Marshall; Gross 1991, Wood et al. 2010. The fundamental cause of the FTC population cycle remains a matter of conjecture. To date, Roland (2005) has provided the clearest evidence that these cycles are driven by a delayed, density-dependent, predator-prey relationship. This relationship varies with local forest structure due to its effects on predator dispersal . Overlaid on this foundational reciprocal feedback relationship are a set of topo-climatic effects involving spring and winter temperatures operating at large and small scales (Daniel and Myers 1995, Roland et al. 1998, Cooke and Roland 2000. In central Alberta, extreme winter temperature fluctuations often result in spatially variable egg mortality , which seems to prevent the emergence of spatially synchronous predator-prey cycling (Hildahl andReeks 1960, Cooke andRoland 2018). In contrast, synchronous herbivore population cycling is often observed in the more moderate climatic regimes of eastern Canada (Cooke et al. 2012) and the northern United States (Witter andKulman 1972, Wood et al. 2010).
Based on spatial analysis of 1950-1984 defoliation data in Ontario, Roland (1993) reported the counterintuitive result that forest fragmentation was associated with an increase in the duration of FTC outbreaks. He suggested that each additional 300 m of forest edge per square kilometer of forest would result in an additional year of outbreak, which is a strong assertion, as it is has been established that while two years of outbreak may not result in appreciable tree mortality, a third year or more of severe defoliation can cause extensive mortality (Churchill et al. 1964, Foster 2017. A small increase in forest fragmentation could thus result in a significant increase in the risk of large-scale forest decline. Subsequent studies in other parts of the insect's range have explored the relationship between forest structure and outbreak parameters, and some have challenged the results of Roland (1993). For example, in Alberta, during the outbreak cycle of 1993-2004, it was found that FTC populations in a 420-km 2 study area exhibited the highest amplitude oscillations in areas that were most heavily forested with trembling aspen (Populus tremuloides Michx.) host trees (Roland 2005). Similarly, in Quebec, Charbonneau et al. (2012) showed that, during the last outbreak cycle of 1998-2003, defoliation impacts were most severe where there was a greater proportion of trembling aspen in the landscape. Furthermore, in the northeastern United States, where sugar maple is the dominant host species, Wood et al. (2010) found that defoliation during the 2003-2006 outbreak tended to last longer in areas where there was more sugar maple. The notion that increasing the amount of host will correspond with increasing insect damage seems intuitive, but conflicts with the results of Roland (1993), thereby producing some lingering uncertainty as to the robustness and importance of the putative effect of forest fragmentation on cycling behavior and outbreak duration. Part of the problem here may be the different ways in which outbreak behavior has been characterized in different studies. Consequently, we begin by illustrating the relationship between various outbreak metrics (i.e., duration, intensity), and ❖ www.esajournals.org show how landscape-driven change in one metric may affect all outbreak metrics simultaneously (Appendix S1: Table S1).
We evaluated the relative importance of FTC host spatial structure (abundance, configuration) and climate on FTC outbreak cycling behavior in a large and heterogeneous landscape located in northwestern Ontario and northern Minnesota: the Border Lakes Landscape (BLL). Prior analyses of the BLL indicate divergent forest landscape legacies across its three principal management zones: older and more conifer-dominated forests within the Wilderness, younger and early forests in an earlier state of succession with more hardwood-dominated forests in managed zones of Minnesota and Ontario, and greater fine-scaled diversity, patch sizes, and edge density within Minnesota relative to Ontario (James et al. 2011, Sturtevant et al. 2014, Robert et al. 2018. It has previously been demonstrated that these differing landscape legacies affect forest insect outbreak dynamics. Robert et al. (2018) showed that spruce budworm (Choristoneura fumiferana, Clemens) outbreaks in the BLL were more intense, more cyclic, and more synchronous in the center of the study area where host spruce and fir trees were more abundant, especially by comparison to the finely fragmented managed forests of Minnesota. This compelling study motivated us to explore whether the spatial legacies of forest management would also be detectable in the demographic histories of other species of outbreaking forest insects, such as the FTC.
Although both the FTC and SBW are lepidopteran defoliators of northern forests, they differ importantly in terms of their preferred hosts. Whereas SBW attacks late-successional coniferous species (i.e., fir and spruce), the FTC is more of a generalist, feeding on hardwood species, especially aspen and sugar maple. Because of these different host preferences, we predicted that the FTC would exhibit the opposite spatial pattern of cycle synchrony from the SBW. Specifically, we expected to identify (1) a breakdown in synchronous cycling within the older and more conifer-dominated forests of the study area core -the Border Lakes Wilderness area; and (2) more synchronous, high-amplitude cycling around the periphery of this core area where abundance of aspen hosts is greater (Appendix S1: Fig. S1).
Further, we predicted that the contrast in the scale of fragmentation provided by Ontario vs. Minnesota management legacies would enable us to confront the hypothesis posed by Roland (1993) within a single ecoregion while minimizing the potentially confounding effect of climate. By studying five full outbreak cycles-not just one partial cycle (as in Roland 2005, Wood et al. 2010, and Charbonneau et al. 2012, or three partial cycles (as in Roland 1993)-and by treating outbreak characteristics as a multivariate response rather than a univariate response, we hoped to explain some of the divergent results being obtained in numerous shorter-term studies. Our overarching hypothesis was that multivariate changes in forest landscape structure affect all aspects of outbreak cycling, because these metrics are functionally dependent on one another (Appendix S1: Fig. S2).

Study area
The BLL contains a large (~20,000 km 2 ) ecoregion that straddles the border between Minnesota (USA) and Ontario (Canada; Fig. 1; Sturtevant et al. 2014, Robert et al. 2018 at the transition between the Laurentian mixedwood and boreal forest regions. The full extent of the BLL includes a 50-km buffer to address potential edge effects related to insect movement (Anderson and Sturtevant 2011). Forest composition is mixed near boreal (Heinselman 1973 Michx.] and trembling aspen), as well as several species near the northern limit of their range, such as white pine (Pinus strobus L.), red pine (Pinus resinosa Ait.), and red maple (Acer rubrum L.).
As described by Robert et al. (2012Robert et al. ( , 2018, the entire BLL shares a common early logging history. Region-wide logging started at the end of the nineteenth century and focused on selective harvest of the big pines (P. strobus, P. resinosa). Logging patterns diverged sharply between Minnesota and Ontario approximately 70 yr ago, coinciding with the advent of mechanized harvesting on both sides of the border. Cut-block sizes in managed Canadian forests are an order of magnitude larger in size than managed American forests, though the rates of harvest have been historically similar (James et al. 2011, Sturtevant et al. 2014. Between these managed regions lies an approximately 10,000 km 2 wilderness area that includes Quetico Provincial Park in Ontario and the Boundary Waters Canoe Area Wilderness in Minnesota, where no timber harvest has occurred since the early 1970s (Heinselman 1996), henceforth referred to as "Wilderness." While forest fires have decreased over the last century (Ahlgren and Ahlgren 2001), enhancing succession toward spruce and fir, periodic defoliation by spruce budworm often kills balsam fir in this region (Frelich and Reich 1995). Actively managed regions of the landscape contain similar forest types as the Wilderness, but have a higher dominance of younger forests, early successional stages (including aspen), due to forestry operations (Wolter and White 2002, James et al. 2011, Robert et al. 2018. Minnesota managed forests have greater local diversity and smaller openings relative to managed forests of Ontario due to the dispersed and small cut-block patterns that occurred there. The divergent landscape structure resulting from distinctly different management legacies is described in more detail in a series of related papers (James et al. 2011, Fig. 1. The Border Lakes Landscape study area located at the border between Ontario (Canada) and Minnesota (United States). Points represent sampling sites for dendrochronological reconstructions of outbreaks that are grouped by the subarea unit of analysis. Subarea labels: First letter is management zone (Ontario, Wilderness, Minnesota); second letter is longitudinal area (East, Central, West); and third letter is latitudinal subarea (North, Central, South). Wolter et al. 2012, Sturtevant et al. 2014, Robert et al. 2018).

Sampling design
We reconstructed FTC outbreak histories using tree-ring analysis of aspen tree cores sampled from 57 sites distributed throughout the BLL (Fig. 1), applying the same site selection criteria and overall sampling design of Robert et al. (2018). Suitable stands were selected when 5-15 host canopy aspen trees were located on mesic flat or mid-slope topographic positions. Two cores were taken from each aspen tree at a height of 1 m. Between 5 and 10 trees were sampled per site (Appendix S2: Table S1).
Groups of three or more neighboring sites were then aggregated into subareas to provide sufficient sample size (n > 15 trees) for outbreak reconstruction (Fig. 1). In total, we created 15 subareas that represent a stratified sampling of the BLL based on forest management zone and geolocation (Appendix S2: Table S2). The goal of this stratification was to account for climatic gradients across the study area. The minimum distance between subarea centroids was 15 km, and the maximum distance was 250 km. Cores were stored in plastic straws and were later mounted, sanded (80, 150, 220, 400 grit size), and stained using Fehling solution (Asherin and Mata 2001).

Outbreak history reconstruction
Forest tent caterpillar outbreaks were reconstructed from trembling aspen tree-ring data by adapting the methods of Cooke and Roland (2007) and Robert et al. (2012Robert et al. ( , 2018. Site-level chronologies were built using tree-ring widths measured using a Velmex UniSlide measuring table with an accuracy of 0.001 mm connected to a computer (Velmex, Bloomfield, New York, USA). Tree chronologies were aggregated into their respective subarea, and we conducted a site-level and subarea-level cross-dating analysis using COFECHA (Holmes et al. 1986) to locate missing or false rings, test for measurement errors, and eliminate problem series. Index chronologies were calculated using the program ARSTAN (Holmes 1992) with a cubic smoothing spline to detrend the series and remove age-related growth trends (Cook and Peters 1981). Spline parameters were set to a 50% frequency response cutoff of 50 yr (Sutton and Tardif 2007). We also used a mean chronology in conjunction with soil moisture index data (Hogg et al. 2013) to assess and partial out the influence of drought on growth (Hogg et al. 2002). From 18 to 42 trees (mean = 25.3 trees; standard deviation [SD] = 6.6 trees) were used to develop mean chronologies for each subarea (Appendix S2: Table S1).
At the level of an individual tree, outbreaks were defined when a tree-ring growth reduction greater than 1.28 SDs from the mean chronology ring width was observed for three or more years (Sutton and Tardif 2007). Outbreak histories, expressed as the percentage of trees classified as affected by an outbreak at a given year and scale, were reconstructed using these parameters and the program OUTBREAK (Holmes and Swetnam 1996). We did not correct the chronologies with nonhost series such as those from jack pine (as in Sutton and Tardif 2007) prior to introducing them to OUTBREAK, because other tree species do not exhibit similar climatic response and because of the high risk of introducing an unwanted semi-periodic signal from other defoliators such as the jack pine budworm, Choristoneura pinus pinus Freeman (OMNRF 2017), or spruce budworm. OUTBREAK parameters used in this study have been used previously to reconstruct outbreak history of FTC in Ontario (Cooke and Roland 2007). The percentage of trees affected by outbreaks each year, at either the site or subarea scales, was used as the response variable for cluster analysis, redundancy analysis (RDA), and analysis of spatial synchrony (see Statistical Analyses).

Predictor variables
Forest landscape structure: general. -Forest landscape structure was quantified using two distinct geospatial datasets derived from multi-temporal Landsat data. The first consists of continuous estimates of basal areas (m 2 /ha) for aspen species and for combined FTC host species (i.e., all nonconifer tree species; Fitzgerald 1995, Wolter et al. 2012. Tree species basal area was quantified at only two points in time (1985 and 2002) due to the greater effort necessary to develop the tree species models (Wolter et al. 2008). The second dataset consists of discrete land-use land cover time-series maps at five-year intervals from 1975 to 2000 (i.e., 6 points in time; Wolter et al. 2012). For both datasets, we did not know a priori the scales at which FTC outbreaks respond to forest conditions. Thus, we quantified our forest landscape structure predictor variables at multiple scales using a range of sample radii (i.e., 0.5, 1, 5, 10, 20, 30 km), where each predictor variable was measured within a circle defined by a given sample radius for each site. Site-level variables were averaged across sites within a subarea (approximately three sites per subareas; Robert et al. 2018). A limitation in any such study is our inability to reconstruct spatial forest histories back to the earliest dates of the tree-ring chronology. The mismatch in dates between our earliest tree-ring data and our earliest satellite data could limit our ability to detect the effects we were seeking to quantify.
To examine the independent effects of FTC habitat abundance and spatial configuration on outbreak dynamics, we divided forest landscape structure into two separate classes of variables, referred to as forest composition variables and forest configuration variables. Composition was quantified using continuous tree species basal area maps, whereas configuration was quantified using discrete land cover maps (Table 1).
Forest composition.-Forest tent caterpillar host forest composition was quantified as the relative basal area of aspen (consistent with Roland 1993 andCharbonneau et al. 2012) and all hardwood species (Nicol et al. 1997), estimated for each sampling radius in both 1985 and 2005. Patterns of tree composition in this study area are also strongly influenced by abundant lakes and wetlands (James et al. 2011). We therefore calculated relative basal area as a percentage of the forested area, and then included percent forest as its own variable, estimated for each of the land cover dates.
Forest configuration. -Forest configuration was quantified using landscape pattern metrics that include edge density of forest/nonforest (ED, km/ km 2 ; Roland 1993) as well as the area-weighted mean patch size of forest openings (AWMPS; McGarigal et al. 2012), and the Shannon-Weaver diversity index (H') of forest classes (Shannon and Weaver 1949) applied to each land cover time interval (Table 1). Each configuration metric was calculated for each of the six available land cover dates (i.e., 1975-2000 at 5-yr intervals) but then averaged across dates in order to capture mean forest structure characteristics of the time period. ED, AWMPS, and H' were calculated using FRAGSTATS (McGarigal et al. 2012).
Climate.-The BLL includes a subtle latitudinal gradient in minimum winter temperature (MWT) such that sites closer to Lake Superior tend to experience a slightly higher MWT (Sturtevant et al. 2014). This gradient is potentially important to FTC because extreme cold or variable temperatures during overwintering can result in outbreak collapse (Blais et al. 1955, Cooke and. Cool springs are also associated with shorter duration outbreaks (Cooke and Roland 2000).
Relevant climate variables were calculated using BioSIM (R egni ere et al. 2014) and included spring degree-days, number of days of spring frost, and winter minimum temperature estimated for each year from 1970 to 2010. Climate data were summarized at the subarea level using the mean of the time series for each predictor variable, and therefore represent the effect of climate normals on spatial dynamics, rather than that of specific weather patterns on temporal dynamics.

Statistical analyses
Time-series clustering.-We used complete linkage cluster analysis (Keogh and Kasetty 2003) applied to time series of raw percentage affected trees to identify subareas with similar outbreak histories. All time series were truncated to the year 1928 as this was the earliest point in time for which we had information on at least five trees per subzone. The resulting set of time series thus span 1928-2006. Clustering was based on the Euclidean distance between time series (Keogh and Kasetty 2003) and used the hclust function in R. Statistical significance of identified clusters was assessed using permutational multivariate analysis of variance using distance matrices as implemented using the adonis function in the vegan (Oksanen et al. 2016) package in R. Here, the matrix of Euclidean distances was used as the response variable and cluster identification was used as a predictor. To determine the degree of periodicity and the primary frequency of oscillation in outbreak chronologies, we conducted spectral analyses on the cluster-level time series (function spec.pgram in R; Daniell kernel; smoothing parameter = 3).
Multivariate analysis.-We adapted the methods of Robert et al. (2018) to model how forest landscape structure and climate affected variation in subarea outbreak chronologies using constrained ordination (i.e., RDA; Legendre and Legendre 2012). In a first application, our multivariate response matrix consisted of the annual percentage of trees affected by FTC outbreak, 1928FTC outbreak, -2006, for each subarea (r = 15 rows; Fig. 1). Potential predictor variables (total of n = 101; Table 1) were of three types: (1) forest configuration (n = 45); (2) forest composition (n = 48); and (3) climate (n = 3). The greater number of predictors in the configuration and composition categories is due to the multiple spatial scales at which these predictors were calculated. Prior to analysis, all predictors were standardized to zero mean and unit variance. We first identified a parsimonious subset of predictors within each predictor type using forward selection and adjusted R 2 as the selection criterion (Dray 2007, Oksanen et al. 2016, Peres-Neto et al. 2006. Forward selection was implemented using the forward.sel function in the packfor package in R. Once predictor subsets were identified for each type, we combined them into a final ordination model. We further reduced the retained predictors by removing those that were identified as collinear based on the variance inflation factor (VIF). Collinearity among predictors can result in biased models and can confound coefficient direction and interpretation. Predictors with a VIF greater than 3.2 were removed from the model (Zuur et al. 2009). Variance inflation factors were calculated using the vif.cca function in the vegan package in R (Borcard et al. 2011).
In a second application of this analytical approach, outbreak regime properties were used as the response matrix in RDA, using the same predictor variables, and variable selection methods outlined above. Outbreak regimes within subareas were characterized using these five attributes: intensity, period, variance, distribution, and synchrony (Appendix S1: Table S1). We predicted that these time-series attributes would vary systematically according to spatial variations in forest landscape structure. Specifically, we predicted that increased abundance and connectivity of host species would result in more frequent and more intense outbreaks with higher amplitude oscillations that resulted in higher variances, less skewed distributions, and higher levels of subarea-scale synchrony.
We chose to analyze both whole time series and outbreak regime attributes to search for general outbreak-landscape pattern associations as well as to test a specific hypothesis about the relationship between outbreak attributes and their mutual dependence on landscape structure (Appendix S1). We predicted that sites would be differentiated according to land-use policies in each of the three management zones (Ontario, Wilderness, Minnesota) and that outbreaks would correspondingly be more intense and better synchronized in areas where host trees were most abundant. Spatial synchrony.-We characterized the spatial synchrony of site-level outbreak time series (n = 57) across the longest interval common to all sites, that is, 1954-2006, using spatial nonparametric covariance functions (SNCFs), as implemented in the ncf package in R (Bjørnstad 2009). Spatial nonparametric covariance function measures the covariance between pairs of sites at varying distances of separation. This analysis typically produces smooth decay patterns in covariance as a function of distance between site pairs. Spatial nonparametric covariance functions were computed for all pairs of sites within each management zone and for all sites within the study area (i.e., a global SNCF). Following Robert et al. (2018), we used site-level data instead of subareas despite relatively fewer replicates per site (5-15 trees) to maximize the spatial resolution of SNCFs within management zones. We predicted that, for any given distance class, there would be lower synchrony (i.e., covariance) within the less-disturbed, less host abundant, Wilderness area than in the managed forest landscapes (Minnesota, Ontario) where aspen host is more abundant.

Outbreak cycling
Between 1876 and 2006, we identified nine reconstructed outbreaks of FTC on trembling aspen where more than 10% of trees exhibited a significant growth reduction within the global average chronology (Fig. 2). These growth reductions occurred in a highly periodic manner, at a frequency of one cycle every 16 yr. The first three outbreak cycles corresponded roughly with outbreaks first identified by Baird (1917) and later confirmed by Ghent (1952) through treering analysis in northwestern Ontario. The cycle peak years were as follows: 1879, 1898, and 1914. The number of tree-ring samples available through the early period 1875-1927 was less than 100, so these data were not considered further in spatial analysis.
The later cycles (i.e., > = 1928), which are the focus of this paper, corresponded with outbreak cycles identified in small-scale tree-ring analyses in Minnesota ) and northwestern Ontario (Cooke and Roland 2007). During the first and last outbreaks of the interval 1928-2006, more than 30% of trees exhibited a growth reduction during the peak year of impact. In contrast, during the middle years , the percentage of trees exhibiting a significant growth reduction never surpassed more than 20% of samples in any one year, even during peak years. During this period, cycling behavior was less well-defined. The average annual percentage of trees affected by FTC was 14.4% (Appendix S2: Fig. S1).
Plot-level outbreak chronologies through 1928-2006 varied considerably, with one subarea (MCS) exhibiting as many as eight distinct cycles of varying amplitude between 1928 and 2006, and one subarea (WSS) exhibiting just two cycles (both severe; Appendix S2: Fig. S1). In sum, there were a total 3-8 cycles observed in the Minnesota subareas, 4-5 in Ontario, and 2-5 in the Wilderness area. The percentage of trees affected by tent caterpillar tended to follow a gamma distribution (Appendix S2: Fig. S3), with three of the four Wilderness chronologies standing out for having highly skewed gamma distributions with many years of zero impact (Appendix S2: Fig. S4).

Cluster analysis
We found several distinct clusters of subarea chronologies such that 70% of time-series variability is captured using six clusters (Appendix S2: Fig. S2). We found that pairs of subarea plots tended to cluster together, with five pairs forming the first five clusters, and the sixth cluster comprised of a set of five chronologies distributed across the study area (Fig. 3). The Minnesota clusters were so distinct that they were the first to group together, with two clusters (1 and 2) alone devoted to distinct variations in Minnesota. Cluster 3 spanned Minnesota (MWN) and the Wilderness region (WWN). Extensive Cluster 4 was rooted in northwestern Ontario (OWN, OWS), but also included subareas from the center of the study area, in Minnesota (MWS) and the Wilderness area (WWS, WSS). Cluster 6 was unique to the Wilderness area. Cluster boundaries tended to overlap in the center of the study area, with clusters 3, 4, and 6 having some representation in the Wilderness area, and clusters 1, 2, and 5 exclusive to managed areas (Fig. 4).
Between five and eight outbreak cycles occurred in each of the six clusters, with ❖ www.esajournals.org substantial variation in both cycle amplitude and intercluster synchrony across the time period (Fig. 5). In particular, the period 1957-1984 was a major source of pattern variation among clusters. Although the global average chronology indicated two partial cycles, IIIa and IIIb, occurring in this window, there was substantial asynchrony in fluctuations during this period, with cycle IV being split into three pulses (IVa, IVb, and IVc), with each pulse affecting disparate parts of the study area. This behavior contrasts with the largely synchronous outbreaks occurring in cycles I, II, and V. The dominant outbreak frequency across clusters was five full outbreaks, yielding a periodicity of 16 yr per cycle, although outbreaks during cycle III were consistently low in intensity. Cluster 1, in the southernmost part of the study area in Minnesota, exhibited a complex pattern of high-frequency oscillations of variable amplitude that made it difficult to identify distinct cycles, with an alternating pattern of three major and four-to-five minor cycles. The alternating pattern of cycle amplitude resulted in a compound periodicity of 26 and 10 yr per cycle. (Periodicities for all series are in Appendix S2: Table S2.) Cluster chronologies were fairly well-synchronized across the study period 1928-2006, with a mean inter-series correlation of 0.42. What the six spatial clusters shared, regardless of location, was relatively high-intensity outbreaks during cycles I and II (in the 1930s and 1950s) and cycles IV and V (in the 1990s and 2000s), leading to the strong positive inter-series correlation of 0.30 outside the 1957-1984 time frame (i.e., 1928-1956 and 1985-2006). Spatial variation among cluster time series was associated with asynchronous features occurring within the 1957-1984 time frame, where one or more cycles were either skipped, split, or shifted in time, leading to a mean pairwise correlation of 0.00 during the low-intensity cycles IIIa and IIIb and the early part of cycle IV.

Multivariate landscape analysis
Redundancy analysis applied to whole time series explained 36.0% (P = 0.001) of the spatiotemporal variation in the annual intensity of outbreaks, with a mix of five forest composition and configuration predictors being identified as significant (Appendix S2: Fig. S5). Climate had no influence on chronology spatial pattern. Reducing the set of forest landscape predictors using VIFs produced an independent set of predictors that explained 16.0% (P = 0.002) of the variation in outbreak intensity (Fig. 6), with no qualitative difference between the five-variable vs. the three-variable ordinations (Appendix S2: Fig. S5). Redundancy analysis axis 1 captured a combined gradient of host abundance and forest gap size with greater host and more closed forest at the left, and open forests and reduced host availability on the right. Redundancy analysis axis 2 captured a gradient in edge density with more forest edges to the bottom of the bi-plot and fewer edges at the top. The three land management zones (O, W, and M) were somewhat separated on the basis of forest landscape structure and outbreak patterning, with Minnesota on the left side of the bi-plot and Ontario and Wilderness to the right. Ontario was more tightly distributed around the center of the bi-plot, with average values for both forest landscape characteristics and outbreak intensity.
Spatial variation in outbreak dynamics is visible at two temporal scales. First, within outbreak pulses (polygons that envelop years in Fig. 6), there is a fast wave-like dynamic where the outbreak peaked first in one set of subareas before spreading to neighboring subareas. The whole sequence typically lasted 3-6 yr from the start to the end of one wave. This is evidence of rapid eruptive growth and spread at smaller spatial and temporal scales. Second, between successive pulses (Ia->Ib->Ic->IIa->IIb->IIIa->IIIb->IVa->IVb->IVc->V), we observed large jumps across bi-plot landscape quadrants (Fig. 6), where one outbreak pulse occurred in one part of the study area (with its characteristic set of forest attributes) and the next occurred in a more distant part of the study area (with a rather different set of landscape attributes).
With 11 different pulses occurring in total, some more intense than others, this is evidence of asynchronous departures from an otherwise regular pattern of synchronized cycling across five area-wide outbreak cycles of 1928-2006. We suspect some other features of outbreak dynamics, such as spatial motion and inertia of an eruptive process, are causing additional separation between clusters and can appear as short-term associations with landscape variables that are not consistent through time. For example, outbreak pulse Ia (1933)(1934)(1935)(1936) was associated with larger forest openings such as observed in OEN, WES, and WEN, whereas pulse IIIb (1970)(1971)(1972) was most associated with WES, where host basal areas measured at Fig. 3. Dendrogram from cluster analysis. Five distinct, well-synchronized pairs, plus a sixth of less synchronous residuals spanning a larger area. See Fig. 1 for subarea labels. 0.5 km in 2005 were highest, meanwhile pulse IVa (1978IVa ( -1981 had its heaviest impact in MCS, where there was a high edge density. In contrast, pulses Ib (1937Ib ( -1940 and IVc (1984IVc ( -1987 were focused on OES, MEN, and MES. The bi-plot indicates that these are subareas with low edge density. In the three-variable model illustrated in Fig. 6, only 16% of the variation in outbreak activity is explained by landscape-meaning that most of this highly dynamic variation in cyclic-eruptive behavior is explained by something that we did not measure. Redundancy analysis on outbreak regime attributes explained only 13.0% (P = 0.02) of the variation in the outbreak spatiotemporal structure (Fig. 7). Climate had no influence on outbreak structure. The only forest landscape variable with any predictive power was host basal area measured at 0.5 km in 1985. The subareas were generally oriented along a gradient that went from high variance and high synchrony in the upper left of the bi-plot to a more skewed distribution of impact intensities in the lower right quadrant. As before, Ontario subareas tended to cluster around the center of the bi-plot, while the Minnesota and Wilderness subareas were more widely dispersed, indicating a wider array of dynamic behavior in these zones. The Wilderness area was more closely associated with the lowerleft edge of the envelope, while Minnesota was more closely associated with the upper-right edge of the envelope. This gradient of behavior was oriented at 45°to the one significant host landscape variable, which is consistent with the weak, albeit statistically significant association Fig. 4. Map of cluster boundaries over the study area. Five of the six clusters were comprised of pairs of subareas located in close spatial proximity to one another. The sixth was an expansive cluster that included parts of Ontario, the Wilderness area, and Minnesota. See Fig. 1 for subarea labels.
between outbreak attributes and forest landscape attributes. Examining correlations between outbreak regime attributes (Appendix S2: Table S3) helps to clarify that there is a relatively strong intrinsic relationship between each attribute, and this strong relationship is, as a whole, weakly influenced by forest landscape structure.

Synchrony within land management zones
As expected, sites closer together tended to exhibit greater synchrony in outbreak chronologies (SNCF; Fig. 8). In all three forest land management zones (Ontario, Wilderness, and Minnesota), spatial covariance dropped off smoothly as the distance between sites grew from 0 to 80 km. The rate of synchrony decay with increasing distance was rapid for the Wilderness area, and much less so for Ontario and Minnesota. The SNCF for the wilderness zone reached zero by 80 km, whereas for the other two zones, even at distances of 200 km, SNCFs seemed to reach a horizontal asymptote near~0.4. This was consistent with the earlier result that cluster boundaries tended to overlap in the Wilderness area, and it is consistent with our initial working hypothesis that the coniferdominated Wilderness area would be a place where cycle synchrony breaks down.

DISCUSSION
Our study investigated the relative importance of landscape forest composition, configuration, and climate on FTC outbreak dynamics within a broad ecoregion containing highly contrasting forest management legacies. This study Fig. 5. Time-series clusters exhibiting weak synchrony within the period 1957-1984 (mean pairwise correlation = 0.00) and strong synchrony outside this period (mean pairwise correlation = 0.30), 1928-1956and 1985-2006 Fig. 1 for subarea labels. represents the longest time period (79 yr and containing at least five periodic outbreaks) for which this question has been investigated for FTC, within a location intentionally selected to increase the probability of identifying the effect of landscape variation while reducing the confounding effects of climate on defoliation patterns.
We found evidence for the existence of five area-wide cycles occurring regularly every 16 yr over the interval . During the intermediate period 1957-1984, outbreak cycle intensity was lower, and cycling was neither regular nor synchronous. Cycle III was clearly broken into discretely regional component pulses (IIIa peaking in 1966 andIIIb in 1972), but even cycles I, II, and IV were partially broken into 2-3 component pulses. Only cycle V was well-synchronized across the study area, and even it was heterogeneous in intensity. Because of the breakdown in synchrony during the 1960s and 1970s, this led to the impression of spatial variation in cycle frequency among subareas and among clusters. From 1928 onward, the number of outbreak cycles varied among clusters from 5 to 7, highlighting the fact that not all cyclic outbreaks grow to an extent that covers the entire study area; many of these cycles are partial outbreaks at best, and they occur at a frequency higher than one every 16 yr in the hardwood-dominated areas of Minnesota.
As expected, we found no evidence for an influence of climate on outbreak dynamics (Figs. 6, 7) at the scale of this study. Climate does nonetheless become important to FTC outbreaks at larger spatial scales (Daniel andMyers 1995, Cooke andRoland 2000). Our analysis was not designed to examine whether low-intensity fluctuations through the 1960s and 1970s might have been a result of prolonged cool weather during Fig. 7. Redundancy analysis bi-plot for outbreak regime characteristics (R 2 = 0.13, P = 0.02). Subareas are black. (Labels in some cases were nudged slightly for enhanced legibility.) Statistically significant landscape variables are depicted in blue. Predictand attributes are plotted in red. The black lines indicate the bounds of an envelope along a gradient of dynamical behavior. See Fig. 1 for subarea labels and Table 1 for landscape variable abbreviations. these decades, but it is a consideration. It is noteworthy, because 1957-1984 is the period where asynchronous patterning drove the clustering pattern, and as Cooke et al. (2012) observed: Low-amplitude fluctuations of weak periodicity tend to be associated with periods of low synchrony, which is consistent with the theoretical predictions of Barbour (1990), and with our own hypothesis.
Our main prediction was that cycle intensity and synchrony would be correlated with one another, and in turn, and would also be correlated with forest species composition and configuration, with more regular outbreak cycles occurring in areas where aspen host trees were most abundant. We found some support for this hypothesis (Fig. 7), but there was a surprising level of complex patterning evident in Fig. 6 that bore no systematic relation to any one landscape variable (see Appendix S2: Table S3 for correlations between outbreak regime attributes, independent of forest landscape structure). As FTC populations cycled slowly across the study area, individual outbreak pulses would emerge and quickly shift as they grew, peaking in an area with a particular landscape attribute, and failing to span the entire study area. Successive cycles would then reach their peak intensity in a different part of the study area, with a different landscape structure. A clear motion of outbreak pulses through time across divergent landscape structures is evident when one traces the 11 successive pulses of outbreaks mapped in Fig. 6. We take this as suggestive evidence that the Border Lakes study area as a whole is sufficiently diverse in forest landscape structure that the cycling dynamics of tent caterpillars are, at best, quasi-synchronous.
We correctly predicted that the Wilderness area is where outbreaks would be least synchronized (Fig. 8); however, we incorrectly predicted that this would result in lower intensity outbreaks. Rather, it appears that while populations in the surrounding landscape are cycling at moderately high amplitude, those waves of outbreak that enter the Wilderness area behave more eruptively, exhibiting a heightened impact, and leading to a tremendously variable set of intensities that are highly skewed in distribution. When populations in the surrounding landscape cycled at low amplitude, such as through the 1960s and 1970s, little or no defoliation was observed in the Wilderness area. It thus appears that FTC dynamic behavior is neither cyclic nor eruptive, but some hybrid mix of these two dynamics.
We expected that by examining multiple outbreak cycles in a large study area with a long history of FTC outbreaks (Sippell 1962, Churchill et al. 1964, Roland 1993, Cooke et al. 2012, we would find strong relationships between outbreak behavior and forest landscape structure. One of our explicit goals was to try to identify stronger and clearer relationships than have been observed in smaller study areas spanning just one cycle in Alberta (Roland 2005), Quebec (Charbonneau et al. 2012), and the American Northeast (Wood et al. 2010).
Although we were able to explain 16-36% of the variation in outbreak activity as a function of forest landscape structure variables, there appeared to be a systematic patterning of impacts in space time that has no obvious explanation: a complex motion where 11 pulse-like waves of outbreak, each occurring every 7 yr on average, are partially synchronized into five large-scale outbreak cycles, each occurring every 16 yr. The result of this fine-scale asynchronous motion-a motion that is far too fast to be explained by slow changes in forest landscape structure-is a set of seemingly strong landscape associations observed in the short run (7-10 yr) that, ultimately, are not durable in the long run (>70 yr). Landscape associations that seem strong when studying eruptive anomalies (individual outbreak pulses) are potentially unrelated to more robust long-term cycling responses to landscape structure, and are therefore likely to be spurious. This provides an extreme caution to landscape disturbance ecologists who study cycling insect systems. Forest managers interested in short-term pest hazard analyses should also take note. Shortterm impact studies are not necessarily a robust basis for land management policy design.
We suggest that our result-a weak association between landscape structure and cycling behavior-is robust and that one of the reasons that published studies disagree on the nature of the relationship between outbreak duration and forest landscape structure is because of the limited extent of temporal sampling that is typically used when studying periodic defoliators. In a stochastic cycling system, where factors such as weather and subtle variations in the natural enemy fauna can strongly affect the expression of any one predator-prey cycle, a single outbreak cycle is not sufficient to investigate an insect's landscape ecology. A larger time-series sample leads to a result that is perhaps unsatisfactory in its statistical significance but is nevertheless more resistant to the sort of errors in model over-fitting that are likely with under-sampled periodic disturbance regimes.
Past studies investigating this question were further limited by applying univariate methods to a multivariate question. In an ideal world, one would select landscapes with orthogonal gradients in the abundance vs. spatial configuration of critical habitats to distinguish the effects of one factor from another, but these two factors are often associated with each other in real landscapes (Fahrig 2017). Likewise, an ideal analysis would examine orthogonal gradients of outbreak attributes, but these characteristics are also related to each other. Given the realities of real landscapes and real outbreak patterns, univariate analyses of either independent (landscape) variables or dependent (outbreak attributes) can be misleading. We suggest that multivariate methods, such as the RDAs applied here, are more appropriate to address the question.
Data limitations are a major caveat with any such study. In particular, historical retrospectives are frequently limited by the unavailability of historical records of landscape structure. In our case, forest compositional variables measured circa 1985 may have a strong influence on subsequent outbreak behavior but may be a poor indicator of forest structures that would have been influencing cycling behavior through earlier decades. Nonetheless, our dataset represents the best available data to ask these large-scale questions.
Finally, although forest landscape structure influenced cycling behavior at the level of land management zone (roughly 100 km in diameter), it appears there may be additional effects manifest at the scale of the study landscape as a whole. Forest tent caterpillar outbreaks were at best quasi-synchronous in this diverse mixedwood landscape, with a complex pulsing and spreading behavior that represented a departure from synchronous cycling. This behavior was surprising and merits further investigation. A constant question in experimental landscape ecology is whether a study was designed at the right scale to evaluate the hypothesis. One is tempted to recommend even larger-scale studies; however, the larger the extent of the study area, the more other factors, such as climate, will enter into play. Given the complexity of the spatiotemporal behavior and the large number of landscape factors at play, applied simulation modeling may be the best way to understand the multiple ways in which forest landscape structure influences spatiotemporal synchronization dynamics.

CONCLUSIONS
Our study represents a comprehensive analysis of the multiple effects of host concentration and forest patch configuration on the long-term dynamics of tent caterpillar cycling. By choosing an experimental landscape that varied significantly in patterns of land-use and forest fragmentation, we were able to examine the independent effects of forest composition and forest configuration that have been identified in earlier studies of FTC landscape ecology.
Our results share some similarity with, but also diverge from, those of Robert et al. (2018) who studied a different cyclic herbivore species -the spruce budworm-in the same BLL and used similar landscape predictors. The consistency is that, at the level of whole landscapes, wherever host forest is more abundant, outbreak cycles are more strongly synchronized. For the FTC, this occurs most clearly in Minnesota and least clearly in the Wilderness area. For the spruce budworm, the opposite is true. Adding nonhost species (conifers in the Wilderness, in the case of FTC; hardwoods in Minnesota, in the case of spruce budworm) tends to push the system in the opposite direction, toward desynchronization. These results demonstrate support for our original hypothesis that any landscape variables that affect cycling dynamics will affect multiple metrics simultaneously, because of their functional dependence on one another. Moreover, this is evidence that host forest species composition can be used as a silvicultural tool to manage pest impacts at the landscape scale, as described by Miller and Rusnock (1993).
Where results differ is in the strength of host landscape influence on cycling behavior. For the FTC, the temporal pattern of fluctuation is highly periodic and the influence of host forest landscape structure is relatively weak. For the spruce budworm, the effect of host forest landscape structure is much stronger. This difference may be a product of the relative strength of the natural enemy interaction thought to be responsible for the intrinsic pattern of cycling-with faster, more robust cycling of FTC (outbreaks lasting usually 1-2 yr and rarely 3-6 yr), and slower, more destructive cycling of spruce budworm (outbreaks commonly lasting 5+ yr per tree and a decade or longer at the landscape scale). If the FTC system is regulated more from the top-down, by natural enemies, and the spruce budworm system is regulated more from the bottom-up, by host plant availability (which diminishes over the course of an outbreak), then this could lead to a divergence in the strength of the influence of host forest landscape structure on spatial patterning.