DOI: 10.5061/dryad.0gb5mkm6q
title: Defoliator outbreaks track with warming across the Pacific coastal temperate rainforest of North America

description: The biogeography of irruptive insect herbivores is determined by host
 availability and climate conditions. As such, outbreak distributions are
 sensitive to climatic change, especially across large latitudinal
 gradients. Here, we investigate the outbreak distributions of two
 understudied defoliators, hemlock sawfly (Hymenoptera; Neodiprion tsugae)
 and western blackheaded budworm (Lepidoptera; Acleris gloverana), that
 have both recently impacted the greatest land area recorded across the
 Pacific coastal temperate rainforest since the establishment of aerial
 survey programs. We compiled polygon-based estimates of insect damage
 collected by aerial observers, forest inventory, and downscaled climatic
 data to develop gridded estimates of bioclimatic conditions across the
 extent of the Pacific coastal temperate rainforest, including the
 continental United States, British Columbia, and Alaska. We leveraged
 these data to develop ensemble machine learning models with the goal of
 predicting the outbreak distribution of each insect. In this manuscript
 we: (1) describe the historical patterns of defoliator outbreaks, (2)
 identify and describe climatic conditions associated with outbreaks in
 both species, and (3) assess whether historic outbreaks have tracked
 geographic shifts in climate conditions across the region. We demonstrate
 that outbreaks of hemlock sawfly and western blackheaded budworm have been
 observed across the Pacific coastal temperature rainforests of North
 America in each decade since the establishment of the Canadian and United
 States aerial survey programs. The distribution of outbreaks by both
 insects were best explained by host availability, a limited range of
 spring, summer, and winter temperatures, and minimum precipitation.
 Finally, we demonstrate that outbreaks have tracked the poleward shift in
 suitable climate over the last century. This study establishes a baseline
 understanding of the climatic constraints and biogeographic patterns of
 historic sawfly and budworm outbreaks across the Pacific coastal temperate
 rainforest and emphasizes the overarching importance of climate in driving
 the irruptive dynamics of these defoliator species.
description: <strong><em>Data</em></strong>
 <em>Aerial surveys</em> We
 compiled aerial detection survey data describing the extent of tree damage
 from insects and diseases. These data were collected by aerial observers
 in fixed-wing aircrafts by the United States Forest Service
 (<em>i.e</em>., Alaska and the western Continental United
 States (CONUS); United States Detection Survey Data) and the British Columbia (BC) Ministry of Forests (British Columbia Overview Survey Data). We calculated the presence/absence of outbreak populations of hemlock sawfly and western blackheaded budworm by intersecting annual polygon and point detection survey data with a 0.25 km<sup>2</sup> grid. Grid cells that contained any visible defoliation were listed as affected. Aerial detection data were available for different time periods in each region of our analysis (Alaska: 1989-2022, BC: 1960-2022, CONUS: 1997-2022). <em>Forest Structure</em> We compiled forest structure data from raster-based estimates of hemlock (<em>spp.</em>) and spruce (<em>spp.</em>) basal area over Alaska and CONUS (per acre; 240 m raster; published in 2012 [Individual Tree Species Parameter Maps [US]) and polygon-based forest inventory records published in 2020 (Vegetation Resources Inventories [BC]). Forest inventory data from the United States were aggregated by upscaling 240 m raster-based estimates of basal area (ft<sup>2</sup>/acre) to 0.25 km<sup>2</sup> and transforming units to m<sup>2</sup>/ha. Inventory data from BC were aggregated by summing the area-weighted mean (per hectare) of forest inventory polygons that intersected each 0.25 km<sup>2</sup> cell. We used total host basal area (hemlock + spruce; m<sup>2</sup>/ha; log-transformed) for all models. <em>Climate</em> We compiled decadal climate data from raster-based estimates of annual and seasonal climate variables (1901-1910…2011-2020; 800 m raster; ClimateNA; Wang <em>et al.</em>, 2016). We selected 14 candidate climate variables: Hogg’s climate moisture index; precipitation as snow (cm); mean annual temperature (˚C); mean annual precipitation (cm); number of frost-free days; spring, summer, autumn, and winter mean temperature (˚C); spring, summer, autumn, and winter precipitation (cm); and minimum winter temperature (˚C). Decadal climate data were intersected with the 0.25 km<sup>2</sup> grid and averaged per grid cell. Climate variables were selected based on hypothesized relationships and tested for inclusion in the final model (Supporting Information). Ultimately, we selected minimum winter temperature, mean spring temperature, mean summer temperature, total summer precipitation, and total winter precipitation for inclusion. <strong><em>Models</em></strong> <em>Ensemble Model</em> We constructed classification ensemble (<em>i.e.,</em> ‘model stacking’) models with the <em>ranger</em> implementation of random forest and the <em>XGBoost</em> implementation of ‘extreme’ gradient boosting. We tuned model hyperparameters based on a grid approach, although neither individual modeling approach is particularly influenced by hyperparameter selection. Based on tuning, we used Gini index as the split criterion for the random forest and grew 300 trees. We fit our gradient boosting models as probability boosted trees with the following parameters: 300 rounds (number of trees); a max depth of 9; η of 0.25 (step size shrinkage), and λ of 0.3 (L2 regularization term). Training data were selected for each ensemble model by randomly selecting 63% of observed presences per decade of observed sawfly and budworm defoliation and an equal number of absences to construct a balanced model. Predictions from each individual model (the random forest and gradient boosting) were averaged to compute a final prediction. We assessed model performances with confusion matrix metrics based on an independently sampled dataset. Since we report model performance across both the full dataset and an independent representative sample, we do not report out-of-bag performance metrics. Further, we used 1-AUC (area under the receiver-operator curve) dropout loss (Bücker et al. 2022) to assess variable importance, and accumulated local effects (Apley and Zhu 2020) to visualize relative mean predictions. <em>Model Workflow</em> We constructed a series of separate ensemble models to predict the biogeographic distribution of hemlock-spruce forestlands and each insect in this study. First, we developed a simple model to predict the presence/absence of hemlock-spruce forestlands for data visualization purposes. We provide the methods and results of this model in Supporting Information. Second, we developed models to predict the spatial extent of outbreak populations of hemlock sawfly and western blackheaded budworm. Models were fit by selecting 63% of presences per decade (1961-1970 onward) and an equal number of absences. Models were constructed 10 times for each insect species and average predictions and performance metrics are reported. Model selection was performed to select the best predictors and predictors were checked for multicollinearity (Spearman pairwise correlation &lt; 0.7). Model selection details are provided in Supporting Information. Finally, based on our model predictions (decadal predictions from 1901-1910…2011-2022) and the thermal climate envelope based on our model results, we explored how the latitudinal range of areas suitable for western blackheaded budworm outbreaks has shifted over the duration of this study. We decided to focus on budworm outbreaks because we have a better historical record of the cyclical outbreaks in this species. <em>Modeling Decisions and Caveats</em> Analyses were restricted to hemlock-spruce defoliation attributed to either hemlock sawfly or western blackheaded budworm. This decision likely results in an underestimate of defoliation because western blackheaded budworm defoliation of spruce is similar to damage caused by other common defoliators, such as western spruce budworm, and thus may be incorrectly attributed; and there are several instances of unspecified ‘unknown’ or ‘defoliator’ reports in forests containing western or mountain hemlock. Further, the reported extent of sawfly and budworm defoliation across southeast Alaska is likely an underestimate due to the limited aerial survey coverage and frequent low cloud ceilings that obscure possible defoliation at higher elevations (see Supporting Information for more details). <em>Model Confidence</em> Constructing spatiotemporally explicit models that predict the likelihood of insect outbreaks across large geographic extents is challenging, especially across areas with low spatial and unequal temporal coverage. We are most confident in our model predictions for Southeast Alaska, Vancouver Island, and the Olympic mountains in Washington. We are the least confident in our model predictions across the coastal mountains of British Columbia and at higher elevations for several reasons: 1) the coastal mountains in BC are sparsely populated and aerial detection survey coverage is extremely sparse; and 2) forestlands across these areas are among the wettest included in our study area. As a result, we are not as confident in predicted relationships across areas where precipitation is extremely high.
description: <em>Packages-</em>All analyses were conducted in R
 (version 4.2.0). Integral packages included: tidyverse (wrangling;
 (Wickham et al. 2019), terra (geospatial formatting; (Hijmans 2024), sf
 (geospatial formatting; (Pebesma 2018), stars (geospatial formatting;
 (Pebesma and Bivand 2023), and mlr (machine learning framework; (Bischl et
 al. 2016).
description: # Defoliator outbreaks track with warming across the Pacific coastal
 temperate rainforest of North America Authors: Michael Howe1,3, Elizabeth
 E Graham2, and Kellen N Nelson1 1 Juneau Forestry Sciences Laboratory,
 Pacific Northwest Research Station, USDA Forest Service, 11175 Auke Lake
 Way, Juneau, AK 99801  2 Forest Health Protection, State, Private, &
 Tribal Forestry, USDA Forest Service, 11175 Auke Lake Way, Juneau, AK
 99801  3 Oak Ridge Institute for Science and Education, Oak Ridge, TN
 37830 Corresponding author:
 [Michael.Howe@usda.gov](mailto:Michael.Howe@usda.gov).  File list:
 defol_data.tiff Ecography_clean.Rmd defol_df.csv grid.hemlock.gpkg
 akse.gpkg na.map.gpkg File Descriptions: Ecography_clean.Rmd -- Analysis
 script. Contains code for models. defol_data.tiff -- Spatial data for
 analysis. defol_df.csv -- Tabular data for analysis. grid.hemlock --
 Spatial data that corresponds to defol_df na.map -- Spatial data that
 serves as the basemap for Figure 1. akse -- Bounding box for plotting
 Southeast Alaska insets. Metadata for defol_data. Data were originally
 processed as 500 x 500 m gridded polygon data. Uploaded data are presented
 in raster form for file size considerations. All raster layers after the
 first four repeat for each decade of analysis: 1901-1910...2011-2020 and
 are summarized below. i5_i5: spatial grid identifier y_y: latitude
 hemlock.spruce_hemlock.spruce: total hemlock spruce basal area (log[m2/ha
 + 1]) hemlock.suitability_hemlock.suitability: predicted hemlock
 suitability (classification model output) Repeating Columns: PPT.sm (mean
 summer precipitation in cm) for 1901-1910... 2011-2020 PPT.wt (mean winter
 precipitation in cm) for 1901-1910...2011-2020 Tave.sm (mean summer
 temperature in ˚C) for 1901-1910...2011-2020 Tave.sp (mean spring
 temperature in ˚C) for 1901-1910...2011-2020 Tmin.wt (minimum winter
 temperature in ˚C) for 1901-1910...2011-2020 bhbw (western blackheaded
 budworm presence/absence) for 1901-1910...2011-2020 hsf (hemlock sawfly
 presence/absence) for 1901-1910...2011-2020 bhbw.pred (predicted outbreak
 distribution of western blackheaded budworm) for 1901-1910...2011-2020
 hsf.pred (predicted outbreak distribution of hemlock sawfly) for
 1901-1910...2011-2020 Metadata for defol_df: i5: spatial grid identifier
 hsf: presence/absence of hemlock sawfly bhbw: presence/absence of western
 blackheaded budworm decade: decade (1901-1910...2011-2020) total:
 hemlock-spruce basal area (log[m2/ha + 1]) PPT.sm: mean summer precip. (m)
 PPT.wt: mean winter precip. (m) Tave.sm: mean summer temp. (˚C) Tave.sp:
 mean spring temp. (˚C) Tmin.wt: min. winter temp. (˚C) hemlock.pred:
 predicted hemlock suitability y: latitude Metadata for grid.hemlock: i5:
 spatial grid identifier x: longitude y: latitude elevation: elevation (m)
 Note: All files must be within the working directory for the analysis to
 work.