Short-Term Increase in Abundance of Foliage-Gleaning Insectivorous Birds Following Experimental Ice Storms in a Northern Hardwood Forest

Leuenberger, Wendy; Cohen, Jonathan B.; Rustad, Lindsey; Wallin, Kimberly F.; Parry, Dylan

doi:10.3389/ffgc.2020.566376

ORIGINAL RESEARCH article

Front. For. Glob. Change, 26 February 2021

Sec. Forest Disturbance

Volume 3 - 2020 | https://doi.org/10.3389/ffgc.2020.566376

Short-Term Increase in Abundance of Foliage-Gleaning Insectivorous Birds Following Experimental Ice Storms in a Northern Hardwood Forest

$\nWendy Leuenberger &#x;$ Wendy Leuenberger¹^*^†

Jonathan B. Cohen¹

Lindsey Rustad²

Kimberly F. Wallin^3,4^†

Dylan Parry¹

¹Department of Environmental and Forest Biology, State University of New York College of Environmental Science and Forestry, Syracuse, NY, United States
²Northern Research Station, USDA Forest Service, Durham, NH, United States
³Rubenstein School of Environment and Natural Resources, University of Vermont, Burlington, VT, United States
⁴Northern Research Station, USDA Forest Service, Burlington, VT, United States

Large-scale disturbances such as ice storms may increase in frequency and intensity as climate changes. While disturbances are a natural component of forest ecosystems, climatically driven alteration to historical patterns may impart fundamental change to ecosystem function. At Hubbard Brook Experimental Forest, NH, experimental ice storms of varying severity were applied to replicate plots of mature northern hardwoods to quantify their effects on forested ecosystems. We assessed ice storm treatment effects on insectivorous foliage-gleaning birds and evaluated insectivore predation on model caterpillars in the understory vegetation. These birds are charismatic, of conservation concern, and are major predators of caterpillars. In turn, lepidopterans are the dominant herbivores in temperate forests and are integral to ecosystem function. We predicted that avian abundance would increase due to additional structural heterogeneity caused by ice treatments, with a concomitant increase in caterpillar predation. Point counts were used to measure insectivorous bird activity in the ice storm experiment plots and additional control plots before and after treatments. We deployed and retrieved plasticine model caterpillars and estimated predation from characteristic marks to these surrogates. Abundance of foliage-gleaning birds was higher in the ice storm plots and birds responded to treatments as a single diffuse disturbance rather than on an individual plot level. All species except one were observed both before and after the ice treatments. Surprisingly, predation on caterpillar models was unaffected by ice storm treatments but rather was a function of caterpillar density. The increase in avian abundance in the ice storm treatment plots corroborates other studies of bird responses to relatively small-scale disturbances in forests and the limited change in species composition was expected given the plot size. We conclude that ice storms may provide beneficial changes for foliage-gleaning birds in the growing season following the disturbance.

Introduction

Concomitant with warming temperatures, climate change is predicted to increase the frequency, intensity, and timing of extreme weather events (Cook and Seager, 2013; Wuebbles et al., 2014). These alterations to natural cycles are directly linked to increases in anthropogenic greenhouse gas emissions (Lambert and Fyfe, 2006; Cook and Seager, 2013) and will likely be more damaging to ecosystems than gradual temperature changes, at least over the short term (Dale et al., 2001; Arnone et al., 2011).

Ice storms are one such event that can have major environmental, societal, and economic consequences as they can span large areas and affect both natural ecosystems and human infrastructure and commerce. Ice storms occur when rain from warm fronts hits colder air closer to the ground (Irland, 1998), resulting in supercooled precipitation and ice accumulation ranging from a light glaze to 10 cm or more (Bruederle and Stearns, 1985; Irland, 1998). From 1832 to 1961, 14 severe ice storms were recorded in New England (Dupigny-Giroux, 2000), a large region prone to frequent icing events. Climate change may increase ice storm frequency and severity as well as shift the current ice belt northward (Dale et al., 2001; Cheng et al., 2011; Klima and Morgan, 2015; Swaminathan et al., 2018).

Ice storms have heterogeneous effects on forested ecosystems as severity depends on stand age, species composition, elevation, proximity to water, slope, aspect, ground surface temperature, precipitation, and wind (Bruederle and Stearns, 1985; Seischab et al., 1993; Irland, 2000; Rhoads et al., 2002; Degelia et al., 2016). Damage wrought by icing ranges from loss of twigs and fine branches to the destruction of the entire crown of large trees (Rogers, 1923; Smith, 2000; Nielsen et al., 2003) with concomitant change to leaf area, canopy cover, and stand height (Rhoads et al., 2002; Weeks et al., 2009). Severely damaged trees can be near trees with moderate or little damage (Rubin and Manion, 2001). Forest structure can be fundamentally altered after ice storms, with an increase in canopy gaps, snags, and a mosaic of downed or partially downed limbs, bent over trees, and coarse woody debris (Smith, 2000; Rhoads et al., 2002; Nielsen et al., 2003), and increased vegetation in the lower strata (Weeks et al., 2009). These changes in vegetation structure and three-dimensional complexity can affect population dynamics of a wide variety of organisms in both positive and negative directions (Blais et al., 2001; Ryall and Smith, 2005; Zhang et al., 2016).

Two trophically linked and ecologically important groups of organisms in temperate hardwood forests are foliage-gleaning insectivorous birds and the lepidopteran caterpillars they consume (Singer et al., 2012, 2014). Ice storms may profoundly affect both groups, through changes in food quality and predation risk for caterpillars (Hirao et al., 2008; Richards and Coley, 2008; Takafumi et al., 2010; Stoepler and Lill, 2013), and through altered prey abundance and habitat complexity for birds. Structural complexity promotes foraging and nesting opportunities for many songbirds that prey on Lepidoptera (DeGraaf et al., 1998; Doran and Holmes, 2005), particularly after disturbances that create openings in an otherwise closed-canopy forest (Blake and Hoppes, 1986; Greenberg and Lanham, 2001). Caterpillar densities and abundance also affect ecological traits of breeding birds in the eastern U.S. including territory size (Marshall and Cooper, 2004), energy intake (Goodbred and Holmes, 1996; Moorman et al., 2007), population growth rates (Townsend et al., 2016), and reproduction (Rodenhouse and Holmes, 1992; Marshall et al., 2002; Newell et al., 2014). Therefore, changes in caterpillar populations are likely to affect many bird populations, including insectivorous bird species of conservation concern (e.g., many passerine Neotropical migrants).

Insectivorous birds exert important top-down control of insect herbivory in forest ecosystems (Bereczki et al., 2015; Muiruri et al., 2015; Whelan et al., 2015). Avian detections are correlated with increased predation on insects (Bereczki et al., 2015) and increased bird abundances resulted in reduced caterpillar population density (Sanz, 2001). Measurements of both predation intensity and avian detections are higher near edges than in interior forest (Bereczki et al., 2015), thus disturbances caused by ice storms that open canopy gaps may drive increases in both variables.

While ice storms occur relatively frequently, these storms are difficult to study due to their stochastic nature, large variability in intensity, and the spatially heterogenous nature of their effects over landscapes (Millward and Kraft, 2004; Wuebbles et al., 2014; Kayler et al., 2015; Degelia et al., 2016). The implementation of experimental ice storms of realistic intensity over a mature hardwood forest canopy in the White Mountains of New Hampshire, USA provided an uncommon opportunity to quantify interactions between insectivorous birds and their lepidopteran prey after glazing events of varying severities. Evaluating the response of these groups to climatically driven disturbances such as ice storms will help us understand potential ramifications of climate changes and determine if management actions are necessary. Our objective was to investigate changes in avian abundance and caterpillar predation rates in response to experimental ice storms. We predicted that increased ice storm severity would result in higher levels of avian abundance and caterpillar predation in the treatment understory.

Materials and Methods

Study Area and Ice Storm Manipulation

To enhance our understanding of the short and long term effects of ice storm on forests, a project was developed at Hubbard Brook Experimental Forest in New Hampshire to experimentally simulate icing events of various realistic severities using high pressure, mechanical pumps to spray water into a hardwood forest canopy under weather conditions where water supercooled, resulting in ice accumulation (Rustad and Campbell, 2012; Campbell et al., 2020; Rustad et al., 2020). The methods are described in Rustad et al. (2020) and Campbell et al. (2020) and summarized in the following section. Ice storm experiment (ISE) plots encompassed ~100 yr old mixed hardwood stands co-dominated by American beech (Fagus grandifolia), sugar maple (Acer saccharum), and yellow birch (Betula alleghaniensis). Ten fixed area plots (20 × 30 m) were assigned to one of five treatments, with two plots in each treatment (Figure 1). Treatments included low (6 mm radial ice accretion), mid (13 mm), high (19 mm), mid×2 (13 mm in two consecutive years), and control (0 mm). Plots encompassed a 5-m buffer around eight 5 m × 5 m subplots and were situated such that their edges were a minimum of 10 m apart. There was some tree, sapling, and shrub removal beginning in July 2015 to create paths on two sides of each plot for vehicles during ice treatments. Ice accumulated over the course of 1–4 h during treatments and remained for five to 13 days, depending on ice thickness and ambient weather conditions. Icing was applied between 18 January 2016 and 11 February 2016. The plots slated to receive treatment in two consecutive years were iced again on 14 January 2017. Data for this study were not gathered in 2017, so the mid×2 were analyzed as replicate mid plots. Icing caused branches to bend, break, or tear from the main stem of the tree and entire trees bent and snapped in the high treatments (Campbell et al., 2020; Rustad et al., 2020). Icing increased deposition of fine and coarse woody debris in the treatment plots (Rustad et al., 2020) and reorganized the canopy with more vegetation at lower strata and a decline in leaf area index following the treatments (Fahey et al., 2020).

FIGURE 1

Figure 1. LiDAR map of ice storm experiment plots with colors representing the target amount of experimentally applied ice accretion. Plots are 20 × 30 m. Hubbard Brook is in the northern part of the figure and the main road is directly south of the plots. Ice storm treatments are located within plots 1–10 and icing occurred in January and February of 2016. Additional control plots include ISEa and plots A100-C250. Circles represent the 15 m point count radii. We visited each point count location six times in 2015 and 2016. The plots were located within Hubbard Brook Experimental Forest in New Hampshire.

In addition to the ten ice storm plots, we sampled six supplementary untreated plots (n = 16 plots total) to address potential edge effects. The edges of untreated plots were located 100 m away (n = 3) and 250 m away (n = 3) from the edge of the nearest ice storm plot (Figure 1). Plots contained comparable vegetation species and structure as the ten ice storm plots. For point counts and canopy cover surveys, we also sampled an additional untreated location between plots 1 and 2 (plot ISEa, Figure 1), for a total of 17 plots as these plots were spaced further apart. Unless otherwise stated, all sampling occurred in both 2015 (baseline) and 2016 (post-manipulation).

Canopy Cover

Canopy cover surveys were conducted in late July to early August. We sampled canopy cover using 15 m transects which began at the point count location (plot center) and extended in each cardinal direction. Canopy cover was sampled at points located every 3 m along each transect, including the center, for a total of 21 points (James and Shugart, 1970). At each point, we used an ocular tube to record presence or absence of canopy cover directly overhead (James and Shugart, 1970). These surveys characterized vegetation damage caused by the ice storms and were used as a covariate for avian detection rates.

Avian Point Counts

We conducted point count surveys to estimate an index of insectivorous bird activity in plots (Petit et al., 1995; Farnsworth et al., 2005). Each plot was visited twice per week from 1 June to 21 June 2016 for a total of six visits per plot per year. All plots were surveyed with 15-m fixed radius point counts. Point counts were conducted from the plot center. We approached the point with minimal disturbance and waited 1 min before beginning a 10-min point count. We identified each individual bird seen or heard during the point count. Multiple individuals of the same species were only recorded when clear evidence of counter-calling occurred (hearing 2+ individuals from different directions in a short span of time). Each round of point counts took place over two mornings, with the 11 ice storm plots on one morning and the six edge effects plots on another morning. All point counts were conducted between sunrise and 08:45 h EST. Point counts were conducted prior to 08:00 h on 96% of counts and in the first hour and half following sunrise (05:05–06:35 h) on 53% of counts. Plot order was reversed on every other visit. Start time, precipitation (none, fog, mist, light rain), cloud cover (0, 25, 50, 75, or 100%), and wind on the Beaufort scale (0–5) were recorded during every count. Counts were only conducted under suitable weather conditions for birds (no heavy rain, winds, or storms). All point counts were conducted by the same observer (WL).

Birds were partitioned into groups based on their foraging behaviors (foliage gleaners, bark foragers, ground foragers, aerial divers) according to descriptions at Cornell University's All About Birds web resource (Cornell University, 2015). We excluded birds that were only detected as flyovers and were not interacting with plots. We retained only foliage-gleaning species in our analyses, as this guild consumes lepidopteran larvae as a major component of their diet (Holmes and Schultz, 1988), alongside other insect groups like Homoptera, Coloeoptera, and Diptera (Robinson and Holmes, 1982). These species generally forage by gleaning, hovering, or hanging on leaves and bark (Robinson and Holmes, 1982).

Caterpillar Predation

To estimate relative predation rates and identify predators of lepidopteran larvae, we used model caterpillars constructed from plasticine (Leuenberger et al., 2019). Plasticine model caterpillars suffer levels of attack similar to real larvae, and thus provide a useful comparison of predation intensity, causation, and impact across treatments (Howe et al., 2009; Bereczki et al., 2014). We made and deployed 50 light green (08 light green Newplast, Newclay Products, Devon, England) geometrid (inchworm) models in each plot (n = 16) for 6 days (Figure 2). There is no consistently-timed peak in caterpillar abundance between late May and late July at Hubbard Brook (Lany et al., 2016), so we deployed caterpillars in late June, which was within the breeding season for all migratory birds at our study area and shortly after peak breeding season for resident black-capped chickadees (Poecile atricapillus) and yellow-bellied sapsuckers (Sphyrapicus varius, Billerman et al., 2020). Accumulated degree days were used to standardize the timing of deployments between years. We obtained daily maximums and minimums for North Woodstock, NH (The Weather Company LLC, 2016). Degree days were calculated beginning on 1 April and included any days with an average temperature ([max + min]/2) above 10°C (Higley et al., 1986). In 2015, we began deployment on 18 June at 262.5 cumulative degree days, and in 2016 on June 23 at 279.7 cumulative degree days.

FIGURE 2

Figure 2. (A) Plasticine caterpillar (left) next to a notodontid caterpillar (right) at Hubbard Brook Experimental Forest, NH, in June 2015. (B) Caterpillar model with evidence of bird predation.

Caterpillar models were 25 mm long with a diameter of 3.5 mm and were glued to leaf petioles or next to leaves in life-like typical day time inchworm poses (Howe et al., 2009). Half of the caterpillars (n = 25) in each plot were deployed on hobblebush (Viburnum lantanoides) and half on American beech. All caterpillars were placed between 0.5 and 2 m above ground, with a minimum of 0.5 m horizontal distance between caterpillars. Caterpillars were spread evenly throughout each plot, typically with 7 or 8 caterpillars per 5 × 5 m subplot depending on the availability of host plants. Some caterpillars were placed in the buffer zone surrounding plots if there were not enough locations within the interior plot to deploy all 50 caterpillars. In 2016, we recorded approximate locations of each caterpillar to the nearest 0.25 m within the plot to allow us to account for spatial autocorrelation (Supplementary Material). Heights were recorded as one of three categories: 0.5–1 m, 1–1.5 m, or 1.5–2 m. Caterpillar locations were discretely marked so as to not provide obvious search cues for predators and checked every other day for evidence of predation attempts (Bereczki et al., 2014).

Any caterpillars with obvious markings from avian predators were removed and placed in 2 mL microcentrifuge tubes for preservation. If the caterpillar was not found on the plant, the ground surrounding the deployment site was searched for 10 min by one person or 5 min by two people. If the caterpillar was still missing after this time, we marked it as “missing.” We continued to briefly check the area on subsequent visits to attempt to recover the caterpillar. All remaining caterpillars were removed after 6 days (Bereczki et al., 2014). Retrieved caterpillars were examined under 1.75× magnification for evidence of predation by birds that could have been missed in the field (Howe et al., 2009).

Data Analyses

Canopy Cover

Canopy cover data were analyzed using logistic regression with a random effect for plot (Bates et al., 2014). Our fixed effects included treatment (Control, Low, Mid, High), site (ISE plots or edge effects plots), and year. We fit models with terms treatment × year, treatment + year, site × year, site + year, each variable as a univariate term, and a null model for a total of 8 models. We considered models with relative model likelihoods of ≥0.125 to have support (Burnham and Anderson, 2002). If multiple models met this criterion, we model-averaged the predicted values using the “AICcmodavg” package (Burnham and Anderson, 2002; Mazerolle, 2016). We simulated residuals using the most complex model (treatment × year) to assess suitability of the distribution and check for overdispersion (Hartig, 2020).

Avian Point Counts

Point count data were analyzed using N-mixture models (Royle, 2004). N-mixture models account for imperfect detection and covariates can be included to model both detection and abundance. Detection is binomially distributed, while abundance is generally Poisson distributed. We could not estimate abundance directly due to the lack of independence among points caused by their small size and tight spacing. Independence among point counts is generally considered to occur when points are located ≥200 to 250 m apart, which diminishes the probability of an individual bird being detected at multiple points (Hutto et al., 1986; Ralph et al., 1993). Therefore, “abundance” as used in this paper refers to a detection-corrected index of the number of individuals using a plot, some of which may have used multiple plots. We fit N-mixture models using the “unmarked” package and R (Fiske and Chandler, 2011; R Core Team, 2019).

We fit the detection component prior to the abundance component of the mixture model. Detection variables included time of day (minutes since midnight), precipitation, sky cover, wind, and canopy cover. All numeric variables were standardized, and detection variables were checked for Pearson's correlations. Variables with correlations >0.7 were not included in the same model (Dormann et al., 2013). Each variable was fit in a univariate model, and in bivariate models with all pairs of these variables, and in an intercept-only null model. When performing model selection for detection, we held the abundance model at one of the most complex possible structures (site × year, Doherty et al., 2012).

The best detection model was used to fit a set of models to assess the effect of treatment, site, and year on avian abundance. Abundance variables included treatment, site (ISE plots or edge effects plots), and year. Models were run for treatment × year, treatment + year, site × year, site + year, each variable as a univariate term, and a null model for a total of 8 models. Each of these sets of abundance variables were run alongside the detection covariate(s) to create the model set. We considered models with relative model likelihoods of ≥0.125 to have support (Burnham and Anderson, 2002). If multiple models met this criterion, we model-averaged the predicted values using the “AICcmodavg” package (Burnham and Anderson, 2002; Mazerolle, 2016). A covariate was considered to have useful predictive value if the 95% confidence intervals for the beta coefficients did not overlap zero. A goodness-of-fit test with 10,000 bootstrapped samples on the site × year abundance model was run to calculate a chi-square test statistic and overdispersion, also known as ĉ (Burnham and Anderson, 2002; Mazerolle, 2016). If ĉ was > 4, we ran the best model with negative binomial and zero-inflated Poisson distributions for the abundance component to see if either method corrected overdispersion (Burnham and Anderson, 2002).

Caterpillar Predation

Caterpillar predation probability was analyzed using logistic regression with a random effect for plot (Howe et al., 2009; Bates et al., 2014). Treatment was included as a covariate, in addition to six categorical and two continuous covariates to account for sources of variability. These covariates pertained to the placement of individual caterpillars (caterpillar height, vegetation), year, and plot location. We used three covariates to split plots according to location, as the density of nearby caterpillars may affect predation: (1) Site: ISE plots and edge effect plots; (2) Location: ISE plots, 100 m edge effect plots, and 250 m edge effect plots; and (3) Periphery: peripheral ISE plots (ISE 1 and 10), interior ISE plots (ISE 2-9), 100 m edge effect plots, and 250 m edge effect plots. Each covariate was run individually, each as an additive variable with treatment, treatment only, and the null model for a total of 14 models. Models with a relative likelihood of ≥0.125 were considered to have support and we model-averaged the predicted models if multiple models qualified as supported models (Burnham and Anderson, 2002). We simulated residuals using the best bivariate model to assess suitability of the distribution and check for overdispersion (Hartig, 2020). Spatial autocorrelation and caterpillar survival over time were calculated (Supplementary Material).

Results

Canopy Cover

Prior to treatment, canopy cover ranged from 71 ± 46% (mean ± SE) to 95 ± 22% for each plot. Following treatment, canopy cover ranged from 52 ± 51% to 100 ± 0% for each plot. The treatment × year model was the only supported model, indicating that the treatments showed different patterns from one another between years (Table 1). Canopy cover decreased with treatment intensity in the second year (Figure 3). Simulated residuals showed that the distribution was a good fit and that there was no overdispersion present in the model.

TABLE 1

Table 1. Model descriptions, number of parameters, and information-theoretic selection criteria for logistic mixed effects models characterizing canopy cover.

FIGURE 3

Figure 3. Canopy cover (model predictions and 95% confidence intervals) before and after ice storm experiment (ISE) treatments. Icing was applied to two 20 × 30 m plots per treatment (n = 4 for medium treatment) during January and February 2016 at Hubbard Brook Experimental Forest in New Hampshire. Additional reference plots were located 100 and 250 m away to investigate edge effects.

Avian Point Counts

We recorded 125 observations of birds from 13 species during point counts (Table 2). Of these, 99 observations from eight species were foliage gleaners with more observations in 2016 (n = 76 observations) than in 2015 (n = 23). All species of foliage gleaners except for the yellow-rumped warbler (Setophaga coronata) were observed in both years (Table 2). One individual could not be identified to species but was foraging in the canopy in a similar manner as foliage gleaners, so was included as a foliage gleaner for analyses. Another unidentifiable bird was foraging on the ground, and thus was excluded from analyses. We recorded between two and 13 observations of foliage gleaners at each point count location over the course of 10 years. We detected a mean of 0.23 ± 0.05 SE detections per count in 2015 and 0.75 ± 0.09 detections per count in 2016. The most common species were red-eyed vireos (Vireo olivaceus; n = 42), black-throated blue warblers (Setophaga caerulescens; n = 16), and black-capped chickadees (n = 10; Table 2).

TABLE 2

Table 2. Total observations of different bird species detected during 15 m radius point count surveys in experimental ice storm plots (n = 11) and controls (n = 6) during 2015 and 2016.

No detection variables had Pearson's correlation values >0.7. The best detection model contained only time. Only seven of the eight models were included in the final model set for abundance as variance could not be estimated in the treatment × year models, likely due to overfitting. Three of the seven models had likelihood values ≥0.125 (Table 3). The abundance components of the three supported models were site × year, site + year, and year (Table 3). Treatment was not included in any of the supported models. The 95% confidence intervals for the beta coefficients did not overlap zero for all coefficients except the abundance intercepts, and year and site in the site × year model (Table 4). Time of day had a positive effect on detection rate in all models (Table 4). Because three models had a relative likelihood ≥0.125, we employed model-averaging on the predicted values to make conclusions. Abundance was highest at the ISE plots in 2016 at 9.14 ± 6.65 SE birds per point, followed by the edge effects plots in 2016 at 4.25 ± 3.12 birds per point, the ISE plots in 2015 at 2.25 ± 1.76 birds per point, and the edge effects plots in 2015 at 2.20 ± 2.11 birds per point (Figure 4). The post-treatment ISE plots differ by <0.2 birds per plot. Time was held constant at the mean value for model-averaging. Confidence limits for all categories overlapped.

TABLE 3

Table 3. Model descriptions, number of parameters, and information-theoretic selection criteria for N-mixture models to model an index of abundance for foliage-gleaning birds.

TABLE 4

Table 4. Regression parameter (Beta) estimates, standard errors (SE), and 95% confidence intervals (CI) for parameters in supported models of an index of avian abundance following experimental ice storms at Hubbard Brook Experimental Forest in New Hampshire during 2015 and 2016.

FIGURE 4

Figure 4. Detection-corrected index of avian abundance with 95% confidence intervals for foliage-gleaning birds in response to ice storm experiment (ISE) treatments. Icing was applied to two 20 × 30 m plots per treatment (n = 4 for medium treatment) during January and February 2016 at Hubbard Brook Experimental Forest in New Hampshire. Additional reference plots were located 100 and 250 m away to investigate edge effects.

We checked goodness-of-fit on the model with time, site × year. Simulated residuals showed that the distribution was suitable and that there was no overdispersion present in the model.

Caterpillar Predation

Predation rates were around 7.5% each year (n₂₀₁₅ = 63 caterpillars, 7.9% of total; n₂₀₁₆ = 59 caterpillars, 7.4% of total). We recovered 792 out of 800 caterpillars per year, for a recovery rate of 99%. Caterpillars not recovered were not found during any checks (n₂₀₁₅ = 5, n₂₀₁₆ = 4), were not found in the final check (n₂₀₁₅ = 3, n₂₀₁₆ = 1), or were accidentally stepped on (n₂₀₁₅ = 0, n₂₀₁₆ = 3). These caterpillars were not included in the analyses. Predation on most caterpillars was identified in the field, with a few examined and confirmed under magnification after the experiment (n₂₀₁₅ = 4, n₂₀₁₆ = 7).

Six models were included as supported models for caterpillar predation by birds (Table 5). These models included site, vegetation type, periphery, location, treatment + height, and treatment + site. Due to this model selection uncertainty, we model-averaged the predicted values to make inferences. No differences among heights were found when all other variables were held constant. More caterpillars were predated on hobblebush than on beech. We then set height to 1 m and vegetation to hobblebush to visually compare the model-averaged predicted values across different locations (Figure 5). Confidence intervals for all groups overlapped. The 100 m and 250 m edge effect plots had a higher predation rate than the ISE plots. The peripheral ISE plots had only slightly more predation than the other eight ISE plots. Predation rates were similar between the 2 years. ĉ of the best model was 1.00, indicating that data were not over dispersed. Caterpillar predation by birds was not autocorrelated at any of the three spatial scales we evaluated (Supplementary Material). Survival rates were higher for caterpillars in the ISE site on all days, without overlap among 95% confidence intervals for day 4 but with overlap for days 2 and 6 (Supplementary Material). Survival rates were lowest in both sites on day 2 and increased on day 4 (Supplementary Material). Survival rates increased for day 6 in the edge effect plots but decreased slightly in the ISE plots (Supplementary Material).

TABLE 5

Table 5. Model descriptions, number of parameters, and information-theoretic selection criteria for logistic regression models of predation probability of plasticine model caterpillars by birds.

FIGURE 5

Figure 5. Model-averaged predation rates for predation by birds on plasticine model caterpillars. There were 50 caterpillars per plot. Plots were located within the ice storm experiment (ISE) plots, some of which were on the periphery (Periph) of the treatments. Other plots were located 100 m or 250 m away from ISE plots to investigate edge effects. Plasticine caterpillars were deployed in late June of 2015 and 2016 and experiments were located within Hubbard Brook Experimental Forest in New Hampshire. Ice storm treatments occurred during January and February 2016.

Discussion

The ice storm experiment resulted in a decreased amount of canopy cover and an increased abundance of foliage-gleaning birds within the ice storm-affected area. Species composition stayed relatively stable, with one additional species observed after the treatments and no other changes. We found some evidence that the increase in abundance was a response to the suite of treatments that created a diffuse disturbance in the roughly 100 × 300 m area containing all eight treated plots and two control plots, although due to model selection uncertainty the evidence is not conclusive. In some regards, these pockets of damaged forest in a matrix of mostly undamaged forest may be akin to windthrow gaps in old-growth or mature forest (Runkle, 1982; Seymour et al., 2002; Campbell et al., 2007) or to a mild natural ice storm or other patchy disturbance, as trees can vary from slightly to severely damaged within a small geographic region (Rubin and Manion, 2001; Faccio, 2003).

Increases in bird abundance following ice treatment corresponds with other regional studies of avian response to small disturbances. Group selection harvests are similar in spatial and temporal scale to these experimental ice storms, as they can be small (0.05–0.80 ha) and the patches can be revisited every 10–20 years (Costello et al., 2000), similar to the recurrence probability of a severe natural ice storm. Group selection cuts of 0.13–0.65 ha in New Hampshire included almost all bird species found in nearby mature forest (Costello et al., 2000). Timber harvests or gaps of 0.5–8 ha often increase avian abundance, species richness, and diversity in forests (Derleth et al., 1989; Forsman et al., 2010). Similar species compositions in control and treated areas are not found in all studies of small disturbances. Neotropical migrants respond differently depending on their habitat associations, as forest-interior species are replaced by interior-edge or edge-open species (Germaine et al., 1997; Moorman and Guynn, 2001; Faccio, 2003; Zhang et al., 2016). We saw little turnover among species, likely because the gap size was smaller than average territory size (Lent and Capen, 1995; Costello et al., 2000). Our study only investigated the first summer following the treatments. In analogous disturbance studies, patterns in avian response to gaps changes with time since harvest with foliage gleaners tending to have a negative relationship with gap age (Forsman et al., 2010).

Tightly spaced gaps and nearby control areas can result in the gaps having a cumulative effect rather than serving as independent replicates, and controls may not be true controls if they are near gaps (Campbell et al., 2007), as exemplified by the control ISE plots that also increased in avian abundance after icing. We saw a small increase in avian abundance at the edge effect plots in 2016. Edge effects from patch cut or group selection harvests have been observed 50–100 m into the forest and effects generally dissipate prior to 200 m (Germaine et al., 1997; Moorman and Guynn, 2001). As such, it is unlikely that the increase we saw is entirely due to edge effects since our plots are up to 250 m away from the ISE disturbance. The increased abundance in 2016 at the edge effect plots may be attributable to environmental stochasticity or to events in localities other than Hubbard Brook Experimental Forest for the migratory species.

We observed a positive relationship between time of day and detection rates in our morning point count surveys. This pattern may be attributable to the end time of our surveys as well as our as well as our species composition. Point counts were completed before 08:00 h so that other researchers could work in the plots without disrupting our point counts. This timing meant that half of our surveys were conducted in the hour and half following sunrise (05:05–06:35 h) when birds are most detectable by point count surveys and did not extend late enough into the morning to observe a decline in detection rates (Robbins, 1981). In addition, we detected more red-eyed vireos than any other species, and this species has a high detection rate throughout the morning (Robbins, 1981). As such, we did not observe the typical decline in detection rate over time seen in other avian point count studies.

We did not see differences in predation of the plasticine model caterpillars. While the density of caterpillars within each ISE or control plot was equivalent, density on a scale relevant to birds was not. Individual plots were at least an order of magnitude smaller than the territory size of the commonly observed species in this study (Steele, 1992; Cimprich et al., 2000; Sillett et al., 2004). Survival rate increased after day 2 in both sites, indicating that birds may quickly learn that model caterpillars are not food and attack fewer as time goes on. Because ISE plots were closely spaced, birds with territories overlapping the ISE plots were likely to encounter a higher density of caterpillars than birds with territories overlapping the edge effect plots. In turn, the 100 m edge effect plots saw lower predation than the 250 m edge effect plots in both years. Some birds may have had territories that overlap both the ice storm plots and the 100 m edge effect plots, but it is less likely that territories also included the 250 m edge effect plots These patterns were backed by our results, as models with categorical terms representing plot density (site, location, periphery) were all supported models and ranked higher than the null or any models including a treatment effect, exemplifying the importance of these plot density terms in explaining predation patterns. The higher survival rates in the ISE plots could be due to plot proximity. Alternative explanations include foraging height, as our model caterpillars were placed at heights of up to 2 m above ground, which was below the subcanopy and canopy levels frequented by the foliage-gleaning species that increased the most between 2015 and 2016 (≥5 more observations; red-eyed vireo, black-capped chickadee, black-throated green warbler [Setophaga virens]; Sturman, 1968, Robinson and Holmes, 1984). Black-throated blue warblers do forage in the shrub level (<2 m. Robinson and Holmes, 1984), and we observed seven in 2015 and nine in 2016. This limited change in black-throated blue numbers may partially explain why we did not see a change in predation rates in their foraging strata. Therefore, observed predation rates may have resulted from an experimental design with tightly spaced plots or the vertical location of our caterpillars in the forest and is likely not representative of foliage-gleaning bird activity. Background density of real caterpillars may also have influenced predation on our model caterpillars, but logistical challenges of enumerating caterpillar density in a mature closed canopy forest are prohibitive so we could not quantify this factor. Similarly, foliage-gleaning birds will also consume spiders and hemipterans (Holmes and Robinson, 1981), which may increase following disturbance and resulting plant growth in gaps (Schowalter et al., 2017).

The foraging ecology of the bird species we observed undoubtedly affected our caterpillar predation study and the response of different species to ice treatments. Red-eyed vireos which forage in the sub-canopy (Holmes and Robinson, 1981) increased in abundance during our study as did black-capped chickadees which are known to search damaged and dead foliage and wood (Robinson and Holmes, 1982). These two species may have increased in abundance due to the altered vegetation structure and increase in dead and damaged wood that resulted from ice storms. Higher predation on caterpillars attached to hobblebush may be a function of an aversion to foraging on beech by several common insectivorous birds at Hubbard Brook (Holmes and Robinson, 1981).

The relatively small plots that formed our experimental units constrained our ability to draw meaningful conclusions about large scale effects of ice storms on avian abundance and foraging. Given the logistical limits of ice storm simulation, it was important to document the response of this feeding guild to the disturbance we created. These spatial challenges are not unique to our study, as Campbell et al. (2007) studied bird abundance and species richness in natural gaps and group-selection timber harvests at Holt Research Forest in Maine and experienced similar mismatches between the experimental design and the appropriate scale to study the organism of interest. In fact, the most frequent natural disturbances in northeastern North America are small canopy gaps typically 0.0024–0.0126 ha (Seymour et al., 2002). These natural gaps are usually caused by tree senescence, wind, mortality from pathogens or insects and have return intervals of 50–200 years (Seymour et al., 2002). Larger, stand-replacing events such as severe ice storms, as well as wind and fire, occur much less frequently (Seymour et al., 2002). As such, our findings provide insight into avian responses to a common size of gap.

Conclusion

Extreme weather events can alter ecosystems but are unpredictable, greatly hampering efforts to understand their effects. Our replicated experimental ice storms provided rare insight into ecological responses after disturbance. Canopy cover decreased and avian abundance increased in the ice storm-affected area following the ice applications. The small size of the treatments meant that we were unable to observe treatment-specific results for birds, but instead showed a general response in the area where the treatments occurred. Few changes in overall species composition were observed. The increase in avian abundance indicates that foliage-gleaning birds may benefit after ice storms, which is likely due to increased structural heterogeneity following such a disturbance. We also showed both the utility and limitations of plasticine model caterpillars in evaluating changes in relative predation intensity.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: Canopy cover data from the Ice Storm Experiment (ISE) plots ver 1. Environmental Data Initiative. https://doi.org/10.6073/pasta/9ba3206dd451ff4547bf50a75813ca12. Plasticine caterpillar predation data from the Ice Storm Experiment (ISE) plots ver 1. Environmental Data Initiative. https://doi.org/10.6073/pasta/da21d491d547876b62d9fee2ed06b635. Avian point count data from the Ice Storm Experiment (ISE) plots ver 1. Environmental Data Initiative. https://doi.org/10.6073/pasta/a3b6394bb75893b076dce82158b74f5c.

Ethics Statement

Ethical review and approval was not required for the animal study because research did not involve an invasive procedure, harm, or materially alter the behavior of an animal under study. No animals were captured or handled.

Author Contributions

WL: conceptualization, data curation, formal analysis, funding acquisition, investigation, methodology, project administration, resources, software, supervision, validation, visualization, writing—original draft, review, and editing. JC: formal analysis, methodology, software, validation, visualization, and writing—review and editing. LR: conceptualization, funding acquisition, project administration, resources, and writing—review and editing. KW: conceptualization, funding acquisition, methodology, project administration, and writing—review and editing. DP: conceptualization, funding acquisition, methodology, project administration, resources, supervision, writing—original draft, review, and editing. All authors: contributed to the article and approved the submitted version.

Funding

Henrietta and John Simeone Fellowship in Forest Entomology (Received by WL; funded WL stipend and research expenses), SUNY-ESF Seed Grant (Received by DP; funded technician stipend), Leroy C. Stegeman Award in Invertebrate Ecology (Received by WL; Funded research expenses), and SUNY-ESF Graduate Student Association Research Grant (Received by WL; funded research expenses).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

Special thanks to S. Hollister, E. Proutey, G. Schwaner, G. Winant, the Ice Storm Experiment team, and the many individuals who assisted with deployment of plasticine caterpillars. Funding for this project was provided by the Edna Bailey Sussman Foundation (WL), Henrietta and John Simeone Fellowship in Forest Entomology (WL), SUNY-ESF Seed Grant (DP), Leroy C. Stegeman Award in Invertebrate Ecology (WL), and SUNY-ESF Graduate Student Association Research Grant (WL). Funding for the main Ice Storm Experiment was provided by NSF (Grant #1457675). Contents of this manuscript also appeared in WL's Master's thesis (Leuenberger, 2017). We appreciate the feedback from the editor and reviewers of this manuscript.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/ffgc.2020.566376/full#supplementary-material

Supplementary Material. It contains methods and results investigating spatial autocorrelation and survival rates of plasticine model caterpillars.

References

Arnone, J. A., Jasoni, R. L., Lucchesi, A. J., Larsen, J. D., Leger, E. A., Sherry, R. A., et al. (2011). A climatically extreme year has large impacts on C4 species in tallgrass prairie ecosystems but only minor effects on species richness and other plant functional groups. J. Ecol. 99, 678–688. doi: 10.1111/j.1365-2745.2011.01813.x