Weather and nest cavity characteristics influence fecundity in mountain chickadees

Nest monitoring

We followed all protocols recommended by the Animal Care Committees of the University of British Columbia and Environment and Climate Change Canada, in accordance with annually renewed permits under the University of British Columbia Animal Care Protocol and Environment and Climate Change Canada’s scientific permit and banding permit number 10365. To find nests, we observed breeding behaviors of adult mountain chickadees and all other cavity-nesting species, and visually inspected tree cavities previously occupied by any species of bird or mammal. We noted which species excavated the nest cavity, and then every 4–7 days we monitored all accessible nests in tree cavities using a video camera cavity monitoring system (for cavities up to 15 m above ground), or a mirror and halogen flashlight from a ladder (up to 6.5 m). We confirmed active nests by the presence of eggs or nestlings. We recorded the date of first egg laid (laying date), clutch size, nest fate (failed or fledged at least one nestling), and brood size (number of young just prior to fledging). We considered clutches complete when we observed the same number of eggs on at least two visits. Nests were deemed successful if one or more chicks were observed to fledge, or were still in the nest at least 17 days after hatching (minimum nestling period for mountain chickadee; McCallum, Grundel & Dahlsten, 2020), with no signs of predation after our last visit with large nestlings. We noted any sign of predation attempts (predator usurped cavity, cavity torn open, or nest material pulled out of the cavity).

We considered second clutches to be nests initiated after either failed attempts, or successful fledging of first broods (Lack, 1954). We were able to confirm a few cases of double brooding (a second brood following a successful first brood, all in 2005 and 2007) by color-banding 115 adults on breeding territories. These adults were captured using mist nets or by covering the nest entrance with a dip net when the adult was inside, banded using a unique combination of three plastic color bands and one numbered aluminium band, and released at the capture site within 10 min. Our field observations and an initial inspection of the data showed a clear bimodal distribution of laying date in 2005 and 2007 (but not in other years). There was a gap in laying around 3 June in both 2005 and 2007, so we considered all clutches laid after 3 June in those years as probable second clutches, and we omitted them from analyses of factors influencing laying date, clutch size, nest survival, and brood size for first nests. In 2002, laying was generally very late (earliest = 21 May, mean = 31 May), and in this year there were 5 clutches laid after 3 June (range 4–26 June) that we classified as probable first clutches.

Nest cavity characteristics

When nesting was finished, we measured cavity height above ground using measuring tape or a 15-m measuring pole, and entrance diameter when it was possible to access the cavity using a ladder. Following the criteria of Edworthy et al. (2018), we classified cavity size as “small” if the entrance measured <3.5 cm in width or if we had observed that the cavity was excavated by a downy woodpecker (Dryobates pubescens), red-breasted nuthatch, black-capped chickadee (Poecile atricapillus), or mountain chickadee (rare, but did occur at 3 cavities that were used for 10 nesting attempts; out of 115 cavities used in 196 nesting attempts with known cavity origin). Cavities were classified as “medium–large” if the entrance measured >3.5 cm in width or if the cavity was excavated by a larger woodpecker (northern flicker; Colaptes auratus, red-naped sapsucker; Sphyrapicus nuchalis, hairy woodpecker; Leuconotopicus villosus, or American three-toed woodpecker; Picoides dorsalis).

Predators and competitors data

To examine site-level factors influencing mountain chickadee fecundity, we established point count stations across our 23 sites (7–32 stations site⁻¹). In continuous forest sites, stations were spaced every 100 m in a grid ≥ 50 m from a grassland or wetland edge (one station ha⁻¹). In forest groves, where it was not possible to establish a grid, we placed point count stations at least 100 m apart. A total of 383 stations were surveyed for mountain chickadees, squirrels, and trees (we collected tree data in 0.04-ha circular plots centered at each point count station; see below); 295 stations were surveyed every year, and 88 had missing data for one or more years. To estimate relative population densities of mountain chickadee and red squirrel, each year from 2000 to 2011 between May and July we conducted 3 rounds of point counts at each point count station. During each 6-minute point count between 0500 and 0930 h, we recorded the species and number of individual chickadees and squirrels seen or heard within a 50-m radius of the station. We detected most individuals by vocalizations, and mapped each individual’s location relative to the centre of the point count station to avoid double counting. To address detection bias by observers across counts, we alternated observers across sites and paired each observer with one of the two principal investigators at least once per season. While it is possible that 6-minute counts could lead to an underestimate of true densities, detection probability was maximised by sampling at point count stations that were typically spaced only 100 m apart. Also, we examined the effect of changes in densities of squirrels and chickadees on fecundity and were therefore concerned only with relative densities and not absolute densities at each site in each year. Further details on point count surveys are provided in Martin & Eadie (1999) and Drever et al. (2008). To generate a measure of squirrel and mountain chickadee density (detections per point count) for each year of the study at each of our 23 sites, we first took the mean number of individuals detected across the three counts at each station, then took the mean across all stations within the site and year. We used a Pearson’s correlation coefficient to examine the relationship between mountain chickadee density from point counts, and the number of nests found in each year (excluding probable second nests).

Prey availability data (site-scale and regional)

Mountain chickadees are opportunistic feeders: in autumn and winter, they live in social groups, feeding on arthropods and conifer seeds (which they cache); in spring and summer, they are socially monogamous, territorial, and feed on insects (including Coleoptera, Lepidoptera, Hymenoptera, and Homoptera), mostly in conifers (McCallum, Grundel & Dahlsten, 2020). Nestlings are fed arthropods (including adults and larvae of mountain pine beetle and western spruce budworm), and at the peak of nestling food demand, adults preferentially take the most readily available prey (Grundel & Dahlsten, 1991; Grundel, 1992; McCallum, Grundel & Dahlsten, 2020). Mountain chickadees exhibit a strong functional response to insect outbreaks, with spruce beetles (Coleoptera: Dendroctonus obesus) comprising approximately 27% of adult diet under outbreak conditions in California, and dozens to hundreds of lodgepole needle miners (Lepidoptera: Coleotechnites milleri) found in chickadee stomach contents under outbreak conditions in Colorado (Telford & Herman, 1963; Baldwin, 1968). In a nest box study in California, they laid larger clutches in years with higher proportions of arthropods (vs. plant material and grit) in the diet prior to the breeding season (Dahlsten & Copper, 1979). Because mountain chickadees forage widely in winter (when they likely accumulate resources to breed, as is the case for a close congener, black-capped chickadee; (Montreuil-Spencer et al., 2019)) but close to the nest in spring (when they begin breeding; McCallum, Grundel & Dahlsten, 2020), we examined prey availability at the regional scale in winter and at the site level in spring. Mountain pine beetle larvae are available to mountain chickadees year-round and were included in both regional and site-level analyses; lepidopteran larvae are only available in spring and summer, and were therefore included only in site-level analyses. To estimate prey availability we used a combination of tree surveys on the ground, and data on the extent of lepidopteran defoliation from the British Columbia Ministry of Forests and Range (Westfall, 2007).

An outbreak of mountain pine beetle occurred across all sites in our study area, with the incidence of attacked pines increasing sharply after 2002; over 95% of the mature lodgepole pine trees (40% of the trees on our sites) were dead by 2005 (Drever, Goheen & Martin, 2009; Edworthy, Drever & Martin, 2011). However, the onset and peak years of the beetle outbreak varied across sites (Drever, Goheen & Martin, 2009; Norris & Martin, 2014). In British Columbia, mountain pine beetles lay their eggs in late summer beneath the bark of living pines, and larvae develop under the bark through the winter and following spring, emerging as adults in summer (Reid, 1962). Outbreaks occur in susceptible forest stands when climatic conditions become favorable for survival of adults and larvae; specifically, a warm August to promote flight by breeding adults followed by a warm winter to enhance survival of larvae, and then low precipitation the next spring (Safranyik, 1978; Safranyik & Carroll, 2006). Beetle larvae available to chickadees during the breeding season were evident the previous summer by the presence of dried resin outflows, or small entry holes (∼2 mm in diameter) in the bark. Each summer, we assessed beetle-infected pine densities for the following spring (chickadee breeding season) by counting the total number of live pine trees ≥ 12.5 cm diameter at breast height (dbh), with bark boring insects, in the 0.04-ha plots centered at each point count station. We subtracted any trees that fell or were cut over the winter, then divided the remaining total by 0.04, and divided again by the number of plots at the site to estimate beetle-infected pines ha⁻¹ at the site level (prey available to mountain chickadees during the breeding season). To assess regional-level mountain pine beetle (prey available to mountain chickadees in each winter), we took the mean number of beetle-infected pines ha⁻¹ across the study area, including only the 295 stations that were surveyed in all years.

Defoliating lepidopterans in our study area include western spruce budworm, which affects living Douglas-fir and spruce, and aspen leaf miner (Phyllocnistis populiella), northern tent caterpillar (Malacosoma californicum pluviale), forest tent caterpillar (Malacosoma disstria), large aspen tortrix (Choristoneura conflictana), and satin moth (Leucoma salicis), which affect aspen (Finck, Humphreys & Hawkins, 1989). Each of these defoliators lay eggs in mid-summer; eggs or early instars overwinter on living leaves or stems/branches of their host tree, and then emerge as larvae in the next spring to summer to feed on the buds and/or new foliage until transforming to adults in mid-late summer (Finck, Humphreys & Hawkins, 1989). Peak outbreak conditions occurred in 2007 for western spruce budworm, and in 2002 for lepidopterans affecting aspen. Under peak conditions, nearly 100% of living host trees at our study area were attacked. To assess lepidopteran density, we combined our field inventory of trees with regional-scale aerial overview survey data for the Cariboo Forest District (excluding Quesnel; BC Ministry of Forests, 1996). We first calculated the density of host trees for each site and year (sum of counts of trees in 0.04-ha circular plots at our point count stations, divided by 0.04, and divided again by the number of plots at the site to give trees ha⁻¹). We used density of living Douglas-fir and spruce trees for western spruce budworm, and living aspen for the remaining defoliators. We then obtained an index of outbreak intensity for each year, from the aerial overview data, by dividing the area defoliated that year by the maximum area defoliated at the peak of the respective outbreak (2002 for aspen defoliators, 2007 for western spruce budworm). We multiplied these indices by the densities of host trees in each site and year to estimate infected trees ha⁻¹ at each site each year for western spruce budworm and aspen defoliators, respectively. Finally, for each site-year combination, we added the estimated density of trees infected by western spruce budworm to the estimated density of trees infected by aspen defoliators, for a total estimated density of lepidopteran-infected trees at each site during the breeding season of each year.

Climate data (regional scale)

Climate measurements were inferred across the study area for each year, by accessing data from the Environment and Climate Change Canada weather station Williams Lake A (WMO ID 71104; 52°10′48″N, 122°03′00″W; elevation 939.7 m; http://climate.weather.gc.ca). We calculated spring temperature as the mean of daily maxima across the period from 7 March to 20 May (the best time window for predicting laying phenology of mountain chickadees, according to Drake & Martin, 2018). Snow depth prior to mountain chickadee breeding was calculated as mean snow depth from 1 January to 30 April. We calculated the number of winter storms (expected to influence laying date and clutch size) as the number of days with rain and/or snow (precipitation) >10 mm between 1 January and 30 April. We calculated the number of spring storms (expected to influence nest survival and brood size) as the number of days with precipitation >10 mm between 1 May and 30 June.

We examined climate variables associated with conifer seed production because mountain chickadees often forage on conifer seeds over winter (McCallum, 1990; McCallum, Grundel & Dahlsten, 2020), and increases in winter food availability can lead to earlier laying dates (e.g., Smith et al., 1980) and increased clutch size in passerine birds (e.g., Broggi et al., 2022). Conifer seed production in dry mixed conifer forests of western North America (in fall) is highly correlated with up to a one-year lag in growing season temperature and soil moisture (Rother, Veblen & Furman, 2015; Petrie et al., 2016); therefore, we calculated 2-year lagged temperature and rainfall for use as predictors in models of mountain chickadee fecundity. We considered rainfall instead of total precipitation because rainfall increases annual soil moisture more consistently than snow due to annual fluctuations in water retention/drainage during snowmelt (Williams, McNamara & Chandler, 2009). We considered 2-year lagged growing season temperature as the mean daily temperature from 1 April through 30 August, 2 years prior to the breeding attempt. Similarly, we calculated 2-year lagged rainfall as the sum of all rainfall from 1 September in year x-3 through 30 August in year x-2, where x corresponds to the year of the breeding attempt.

Statistical analyses

We used mixed-effects models to test a priori hypotheses about the influence of regional, site-level, and nest-level factors on four components of mountain chickadee fecundity (response variables for each nest were: laying date, clutch size, nest survival, and brood size; Table 1 and Table S1). Probable second nests (first egg laid after day 150 in 2005 or after day 160 in 2007) were excluded from all models. We elected to run separate models predicting each component of fecundity (rather than a single path analysis) because many nests had incomplete fecundity data, resulting in different but overlapping datasets for each response variable. We assumed that the nests for which appropriate fecundity data were available were independent and representative of all nests present in each year. At the regional scale, we constructed two models for each fecundity response variable, and considered models to be plausible (supported by the data) if their ΔAICc values were within 2 values of the top model. At site- and nest-scales, we constructed one model for each fecundity response variable. At all scales, fixed effects were considered to have statistical support from the data if the 95% confidence intervals on parameter estimates did not overlap zero, and if they were included in a model with ΔAICc values >2 points below the null (intercept-only) model. We predicted that the number of storms would affect the correlation between annual fecundity and spring temperature such that winter storms would impede the ability of females to lay earlier and larger clutches in warmer springs, and spring storms would interfere with feeding nestlings leading to greater partial nest loss (de Zwaan et al., 2020). Therefore, we included interaction terms between spring temperature and storms for models explaining variation in each fecundity response variable.

All data analyses were conducted in R version 3.5.3 (R Core Team, 2019); packages listed in Table S1). Continuous predictors were scaled to have a mean of 0 and standard deviation of 1. Due to multiple measurements of fecundity within cavities (nests) and sites across years, the error term of nest and site crossed with year would have appropriately reflected the structure of the data. However, including site as a random effect resulted in singularity in some models, likely as a result of overfitting. Therefore, we applied a model selection approach to choose the model with the random effects term that balanced predictive accuracy and overfitting/type I error (Bates et al., 2015; Matuschek et al., 2017). We excluded the random effect of site because where the models converged and did not result in singularity, a comparison between top models for laying date, clutch size and brood size, with random effects of (1) year only and (2) site crossed with year, the model with year only was better supported by the data (ΔAICc <2). We lacked sufficient data to include a random error term to account for within-nest variation (at the nest scale), and therefore were unable to compare models across scales, so we limited our model selection approach to within scales only. For nest survival, some models with a random effect of year resulted in singularity, and removing the random effect resulted in very similar parameter estimates. Thus, we included only year as a random effect in models of laying date, clutch size, and brood size, and no random effects in models of nest survival. We examined the correlation coefficients of predictors using correlation matrices of fixed effects for each model, and found no significant correlations (r < 0.7). We modeled laying date (continuous response variable) using general linear mixed models with a Gaussian error structure, and examined residuals and Q-Q plots to check model fit and assumptions (Pinheiro et al., 2021). We modeled clutch size and brood size (under-dispersed counts), using generalized linear mixed models with COM-Poisson error structure to estimate a dispersion parameter to handle underdispersion (Sellers, Borle & Shmueli, 2012; Brooks et al., 2017). We modeled daily nest survival (binary) using the logistic exposure method described by Shaffer (2004), a variation on logistic regression in which the link function accounts for the number of exposure days between observer visits to the nest. In modeling brood size, we restricted the dataset to nests that were successful. In modeling nest survival and brood size, clutch size was included as a predictor, so data were restricted to nests with complete clutches of eggs observed (difficult to obtain given the configuration of the nest cups within tree cavities and the chickadees’ habit of covering their eggs before leaving the nest). Results are reported as mean ± standard deviation.

Annual variation in fecundity

In 12 years across 23 sites, we studied 513 mountain chickadee nests in tree cavities primarily excavated by woodpeckers and nuthatches. Incubation period was 13 ± 1 days (n = 60 nests) and nestling period was 19 ± 1 days (n = 14 nests). Initiation of egg laying ranged annually over a 10 week period (earliest clutch initiation was on 1 May in 2005 and latest on 9 July in 2007 (Fig. 2A). In most years, laying was finished by 9 June; however, in 2005 and 2007, there were two peaks in laying date—one in mid-May, and the second in late June. At least 32 of the late nesting attempts (and none of the early nesting attempts) were second clutches in a tree or vicinity where mountain chickadees had nested earlier in the season (confirmed by banding to be the same pair in six cases). Overall, mean laying date was 20 May for 129 known and probable first clutches, and 23 June for 32 known and probable second clutches. Mean hatch day was 6 June ± 8 days for 45 known first clutches, and 32 days later, 8 July ± 7 days for 10 known second clutches. On two occasions, parents were observed feeding fledglings of the first brood while attending the second nest.

Clutch size was 6.81 ± 1.25 eggs (range: 3–10) for known and probable first clutches (n = 189; Fig. 2B). Known and probable second nests had a clutch size of 6.24 ± 1.23 eggs (range: 3–10) (n = 33). We determined nest fate for 208 nests, of which 142 (68%) fledged at least one nestling (67% of 184 first nests and 75% of 24 second nests); nest success was highest in 2005 (88% of 33 first nests, 80% of 20 second nests) and lowest in 2011 (29% of 7 nests) and 2003 (40% of 15 nests; Fig. 2C). We found evidence of predation at 45 of the 66 failed nests (nest was empty before 14 days post-hatch and the cavity was usurped by a predator, torn open, or nest material removed). At successful nests, mean brood size of first nests was 5.48 ± 1.55 fledglings (range: 1–9) and varied annually, with mean brood size highest in 2007 (6.30) and 2005 (6.18) and lowest in 2003 (4.50; Fig. 2D). Mean brood size for second nests was 5.21 ± 1.37 fledglings (range: 1–7). The total number of probable first nests found in each year was positively correlated with the annual mean number of mountain chickadee individuals detected per ha in point counts (r = 0.70, df = 10, p = 0.01).

Factors influencing fecundity in first nests

Models at the regional scale explained some of the variation in laying date and nest survival of mountain chickadees; models at site scale (insect availability, conspecifics, predators) performed poorly at explaining variation in fecundity (no better than the null model); and models at nest scale explained some of the variation in clutch size, nest survival, and brood size (Table 2, Fig. 1). Earlier laying was associated with warmer spring temperatures; however, this relationship was reduced with increasing number of storms from Jan 1 through Apr 30 (Tables 2 and 3; Figs. 1 and 3). Clutch size declined with laying date both within and across years (Tables 2 and 3; Figs. 1 and 4). Daily nest survival (DSR) increased with annual temperature 2 years prior to the nesting attempt at the regional scale, and increased with clutch size but declined with cavity size at the nest scale (Tables 2 and 3; Figs. 1 and 5). Extended over the 39-day nesting period from clutch initiation to fledging, probability of nest survival (DSR³⁹) ranged from 21% in a large cavity with a clutch size of 4, to 95% in a small cavity with a clutch size of 9, and was 88% for the mean clutch size of 6.81 in (typical) small cavities. Brood size (in successful nests), was strongly and positively predicted by clutch size only, and year effects were negligible in this relationship (Tables 2 and 3; Figs. 1 and 6).

Three models included statistically significant effects related to climate and predators, but were no better than the null models at explaining variation in fecundity. In the regional climate models, clutch size and nest survival increased with spring temperature except when there were more winter storms (Tables 2 and 3). At the site level, brood size declined with increasing squirrel density (Table 3). None of the chickadee fecundity measures were predicted by our measures of conspecific densities, cavity height, or snow depth in any of the models (Table 3).